{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# For-loops with theano or confrontational therapy for `theano.scan`-phobia\n",
    "\n",
    "In this notebook, we will gently approach the function `theano.scan` from known concepts in python. `theano.scan` is a very general metafunction that permits the evaluation of loops fully within a computational graph. This implies that gradients can be automatically calculated from the resulting expressions. \n",
    "\n",
    "At the end of this tutorial, you should be able to understand and manipulate code like this (taken from the RNN tutorial)\n",
    "\n",
    "```python\n",
    "def step(x_t, h_tm1):\n",
    "    h_t = T.nnet.sigmoid(T.dot(x_t, self.W_in) + \n",
    "                            T.dot(h_tm1, self.W_rec))\n",
    "    y_t = T.nnet.softmax(T.dot(h_t, self.W_out))\n",
    "    return [h_t, y_t]\n",
    "            \n",
    "[h_vals, y_vals], _ = theano.scan(fn=step,        \n",
    "                                  sequences=X.dimshuffle(1,0,2),\n",
    "                                  outputs_info=[h0, None])\n",
    "```\n",
    "\n",
    "\n",
    "`theano.scan` has many arguments and is quite formidable at first sight. We will introduce them one by one via two specific derivations of `theano.scan` which have pure-python counterparts: `map` and `reduce`.\n",
    "\n",
    "We start with a quick refresher of functional programming concepts in pure python. Then we explore the theano equivalents. After this we will be equipped to tackle `theano.scan` in its full generality and at the end we'll code a small dynamical system."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "%matplotlib inline\n",
    "import theano\n",
    "import theano.tensor as T\n",
    "import numpy as np\n",
    "import matplotlib.pyplot as plt\n",
    "\n",
    "theano.config.floatX = 'float32'\n",
    "theano.config.optimizer='fast_compile'"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Micro-reminder of functional programming keywords\n",
    "We begin by recalling what `map` does: It applies a given function to all elements of an iterable sequence. We can use it, for example, to lift a sequence of integers to their squares."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "[0, 1, 4, 9, 16, 25, 36, 49, 64, 81]"
      ]
     },
     "execution_count": 2,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "def square(x):\n",
    "    return x ** 2\n",
    "map(square, np.arange(10))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Relatedly, the function `reduce` accumulates (or 'folds') a sequence by iterative application of a binary function. The endpoint operation is done with a neutral element or the endpoint of the list. For instance, we can add up a sequence of numbers"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "45"
      ]
     },
     "execution_count": 3,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "def add(a, b):\n",
    "    return a + b\n",
    "reduce(add, np.arange(10), 0)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We can also apply it after modifying the sequence with a call to `map`, giving rise to an instance of the widely employed 'map-reduce' framework:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "285"
      ]
     },
     "execution_count": 4,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "reduce(add, map(square, np.arange(10)), 0)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Make a verbose form of `add` called `verbose_add`, which prints its arguments when it is called. Map a function that takes a sequence element and returns a singleton list (`lambda x: [x]`) to a sequence of integers and reduce it using `verbose_add`. Try specifying different neutral elements or leaving them out."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "def verbose_add():\n",
    "    pass\n",
    "\n",
    "#reduce(lambda ...)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "a=[100], b=[0]\n",
      "a=[100, 0], b=[1]\n",
      "a=[100, 0, 1], b=[2]\n",
      "a=[100, 0, 1, 2], b=[3]\n",
      "a=[100, 0, 1, 2, 3], b=[4]\n",
      "a=[100, 0, 1, 2, 3, 4], b=[5]\n",
      "a=[100, 0, 1, 2, 3, 4, 5], b=[6]\n",
      "a=[100, 0, 1, 2, 3, 4, 5, 6], b=[7]\n",
      "a=[100, 0, 1, 2, 3, 4, 5, 6, 7], b=[8]\n",
      "a=[100, 0, 1, 2, 3, 4, 5, 6, 7, 8], b=[9]\n"
     ]
    },
    {
     "data": {
      "text/plain": [
       "[100, 0, 1, 2, 3, 4, 5, 6, 7, 8, 9]"
      ]
     },
     "execution_count": 6,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# %load scan_tutorial_solutions/sol_reduce_verbose_add.py\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Observe that the neutral element is fed at the beginning and that the accumulated values are fed into the left-hand argument. (**This will be the other way round in Theano.**) \n",
    "\n",
    "\n",
    "Note that `map` can actually take multiple sequences as input"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "a=0, b=100\n",
      "a=1, b=101\n",
      "a=2, b=102\n",
      "a=3, b=103\n",
      "a=4, b=104\n",
      "a=5, b=105\n",
      "a=6, b=106\n",
      "a=7, b=107\n",
      "a=8, b=108\n",
      "a=9, b=109\n"
     ]
    },
    {
     "data": {
      "text/plain": [
       "[100, 102, 104, 106, 108, 110, 112, 114, 116, 118]"
      ]
     },
     "execution_count": 6,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "map(verbose_add, np.arange(10), np.arange(100, 110))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Some functions with sequence output\n",
    "We will now work with functions that could be called 'reduce with history', but that is a cumbersome name.\n",
    "\n",
    "In the following you will be writing a lot of functions that are similar to each other, but with increasing complexity. Feel free to copy and paste from your previous work when appropriate.\n",
    "\n",
    "Next up is a function that cumulatively adds all integers up to a certain point, starting from 0"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "[0, 1, 3, 6, 10, 15, 21, 28, 36, 45]"
      ]
     },
     "execution_count": 7,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "def add_ints(n):\n",
    "    output = [0]\n",
    "    for i in range(1, n):\n",
    "        output.append(output[-1] + i)\n",
    "    return output\n",
    "\n",
    "add_ints(10)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Write a function that adds up the first n squares cumulatively and returns all of them."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "def accumulate_squares():\n",
    "    pass"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 101,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "# %load scan_tutorial_solutions/sol_pure_python_accumulate_squares.py\n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "[0, 1, 5, 14, 30, 55, 91, 140, 204, 285]"
      ]
     },
     "execution_count": 12,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "accumulate_squares(10)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Write a function that generates n elements of the Fibonacci sequence, using the recurrent formula $x_{n+2} = x_{n + 1} + x_n$,with $x_1 = 1$ and $x_0 = 0$."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "def fibonacci(n):\n",
    "    pass"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "# %load scan_tutorial_solutions/sol_pure_python_fibonacci.py\n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "[0, 1, 1, 2, 3, 5, 8, 13, 21, 34]"
      ]
     },
     "execution_count": 16,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "fibonacci(10)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now write a version of `accumulate_squares` that takes an arbitrary sequence as input"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "def accumulate_squares():\n",
    "    pass"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 23,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "# %load scan_tutorial_solutions/sol_pure_python_general_acc_squares.py\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 22,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "[0, 4, 20, 56, 120, 220]"
      ]
     },
     "execution_count": 22,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "\n",
    "\n",
    "accumulate_squares(np.arange(0, 12, 2))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Take the next generalization step by providing a function that cumulatively combines any given sequence, similar to `reduce` but keeping the history. Then rewrite `add_squares` and `add_ints` using it. Make sure the left-hand argument of the binary function receives the current sequence value and the right hand argument receives the accumulated value. Also provide an argument `starting_point` from which to initiate the recursion but which will not appear in the output."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "def accumulate():\n",
    "    pass"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 28,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "# %load scan_tutorial_solutions/sol_pure_python_general_accumulate.py\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 27,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "a=0, b=4\n",
      "a=2, b=4\n",
      "a=4, b=6\n",
      "a=6, b=10\n",
      "a=8, b=16\n"
     ]
    },
    {
     "data": {
      "text/plain": [
       "[4, 6, 10, 16, 24]"
      ]
     },
     "execution_count": 27,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "accumulate(verbose_add, np.arange(0, 10, 2), 4)    "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "[0.0, 1.0, 5.0, 14.0, 30.0, 55.0, 91.0, 140.0, 204.0, 285.0]"
      ]
     },
     "execution_count": 12,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "accumulate(lambda a, b: a ** 2 + b, np.arange(10), 0.)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## The theano versions of `map` and `reduce`\n",
    "The theano versions of `map` and `reduce` work very similarly to their python counterparts, except that all entries can be symbolic. Let us first take a look at the docstrings and then try out some examples"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### `theano.map`"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "theano.map?"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Note that the function `map` returns `expression, updates`, of which we right now only need `expression`. The following code creates a theano function that squares the first `n` integers."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 60,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([  0,   1,   4,   9,  16,  25,  36,  49,  64,  81, 100], dtype=int32)"
      ]
     },
     "execution_count": 60,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "n = T.iscalar()\n",
    "map_square_expr, updates = theano.map(square, sequences=T.arange(n))\n",
    "\n",
    "f_map_square = theano.function([n], map_square_expr)\n",
    "\n",
    "f_map_square(11)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Write a function that maps the function `square` to any input array."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "# ...\n",
    "# f_square_arr = theano.function(...)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 32,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "# %load scan_tutorial_solutions/sol_theano_map_square.py\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 31,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([   0.,    4.,   16.,   36.,   64.,  100.,  144.,  196.,  256.,  324.], dtype=float32)"
      ]
     },
     "execution_count": 31,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "f_square_arr(2. * np.arange(10).astype(np.float32))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now, using `theano.map`, write a function that takes two arrays and adds them elementwise. Use the function `add`. Make sure that the function returns what you expect."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "# ...\n",
    "# f_add = theano.function(...)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 36,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "# %load scan_tutorial_solutions/sol_theano_map_add_two_seq.py\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 35,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([ 100.,  102.,  104.,  106.,  108.,  110.,  112.,  114.,  116.,  118.], dtype=float32)"
      ]
     },
     "execution_count": 35,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "\n",
    "\n",
    "f_add(np.arange(10).astype('float32'), np.arange(100, 110).astype('float32'))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Take another look at the docstring of `theano.map` and identify the argument `non_sequences`. This is a list of auxiliary expressions that can be passed to the scan function. While `theano.maps/reduce/scan` sees all expressions outside its arguments, it is preferable to pass any addition expression that will be used via `non_sequences`. We can show this using the example of a polynomial of fixed coefficients whose list of monomials we will generate. \n",
    "\n",
    "A polynomial in $x$ can be written as $p(x) = c_0 x^0 + c_1 x^1 + c_2 x^2 + ... + c_N x^N$. We will keep the coefficients $c_i$ in an array `coefs` and define the symbolic variable `x`. We choose the coefficients to be $c_i = \\left(\\begin{array}{c}4\\\\i\\end{array}\\right)$, such that $p(x) = (x + 1)^4$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 42,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "coefs = T.constant(np.array([1, 4, 6, 4, 1]).astype('float32'))\n",
    "x = T.fscalar()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Using `theano.map` and an appropriate function, create an expression that outputs a list of monomials in the symbolic variable `x` with coefficients from `coefs`."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 39,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "def generate_monomial():\n",
    "    pass\n",
    "\n",
    "#...\n",
    "# f_monomials = theano.function(...)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 45,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "# %load scan_tutorial_solutions/sol_theano_map_univar_polynomial.py\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 44,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "81.0"
      ]
     },
     "execution_count": 44,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "f_monomials(2.).sum()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 50,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "y = T.fscalar()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now do the same for two variables `x, y` in order to calculate the list of monomials for `(x + y) ** 4`"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 53,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "def generate_double_monomial():\n",
    "    pass\n",
    "\n",
    "# f_binomial = theano.function(...)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 54,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "# %load scan_tutorial_solutions/sol_theano_map_binomial.py\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 52,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "81.0"
      ]
     },
     "execution_count": 52,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "\n",
    "f_binomial(1., 2.).sum()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We can also take derivatives through these polynomial expressions. Let's try this for the univariate polynomial first. We have $p'(x) = 4(x + 1)^3$ and so expect, e.g. $p'(1) = 4\\times 2^3 = 32$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 55,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array(32.0)"
      ]
     },
     "execution_count": 55,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "univar_polynomial_gradient = T.grad(monomial_expr.sum(), wrt=x)\n",
    "\n",
    "f_univar_pol_grad = theano.function([x], univar_polynomial_gradient)\n",
    "\n",
    "f_univar_pol_grad(1.)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now do the same for the multivariate polynomial. We expect \n",
    "$p(x, y) = \\left(\\begin{array}{c}4(x + y)^3\\\\4(x + y)^3\\end{array}\\right)$\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 58,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "[array(4.0), array(4.0)]"
      ]
     },
     "execution_count": 58,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "multivar_polynomial_gradient = T.grad(double_monomial_expr.sum(), wrt=[x, y])\n",
    "\n",
    "f_univar_pol_grad = theano.function([x, y], multivar_polynomial_gradient)\n",
    "\n",
    "f_univar_pol_grad(0.5, 0.5)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### `theano.reduce`"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now let's take a look at the docstring of `theano.reduce`."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 23,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "theano.reduce?"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "As you can see, it takes 3 compulsory arguments - the function `fn`, a list of `sequences` and `outputs_info`. It will become extremely important to understand the role of the latter. It has several functions, the first of which we shall concentrate on now: Providing initial values (or the 'neutral element' for the reduce operation).\n",
    "\n",
    "We begin by using `theano.reduce` to sum a sequence of integers."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 61,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array(45, dtype=int32)"
      ]
     },
     "execution_count": 61,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "int_sum_expr, updates = theano.reduce(add, sequences=T.arange(n), outputs_info=np.int32(0))\n",
    "\n",
    "f_int_sum = theano.function([n], int_sum_expr)\n",
    "f_int_sum(10)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Note that it is highly important to get the data types right - this implies being aware of potential implicit type casts within functions.\n",
    "\n",
    "Now write a function that adds all squares of the first `n` integers."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 65,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "def add_square_a_to_b():\n",
    "    pass\n",
    "\n",
    "# f_square_sum = theano.function ..."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 104,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "# %load scan_tutorial_solutions/sol_theano_reduce_sum_square_ints.py\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 64,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array(285, dtype=int32)"
      ]
     },
     "execution_count": 64,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "\n",
    "\n",
    "f_square_sum(10)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now we would like to reduce a function with multiple outputs. Find out how to employ reduce to return the sum of the first `n` integers *as well as* the first `n` squares"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 45,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "[array(45, dtype=int32), array(285, dtype=int32)]"
      ]
     },
     "execution_count": 45,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "def add_squares_and_ints(i, cur_int_sum, cur_square_sum):   # to be removed\n",
    "    return i + cur_int_sum, i ** 2 + cur_square_sum\n",
    "\n",
    "both_expressions, updates = theano.reduce(add_squares_and_ints,\n",
    "                                          sequences=T.arange(n),\n",
    "                                          outputs_info=[np.cast['int32'](0)] * 2)\n",
    "f_both = theano.function([n], both_expressions)\n",
    "\n",
    "f_both(10)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now look at the following function"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 68,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "def return_args_and_accumulate_squares_and_ints(i, cur_int_sum, cur_square_sum):\n",
    "    return i, cur_int_sum, cur_square_sum, i + cur_int_sum, i ** 2 + cur_square_sum"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Use this function in `theano.reduce` to obtain exactly the same result as before. Note that the first three outputs of the function should be considered as useless and should be ignored. To do this, you need to modify `outputs_info` and place `None` at every place where you are not interested in feeding the output of the function back.\n",
    "\n",
    "We have just discovered the second meaning/use of `outputs_info`: Aside from containing initial values, it also regulates which values of the function are fed back. Thinking about this a little, the two concepts go hand in hand: Function outputs that are fed back into function require arguments to take those values. These arguments need to be fed with initial values at the beginning of the loop.\n",
    "\n",
    "Note that the final output will also contain the final outputs of the function that were not fed back."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "# ...\n",
    "\n",
    "#f_both = theano.function ..."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 73,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "# %load scan_tutorial_solutions/sol_theano_reduce_multi_output.py\n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 72,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "[array(9, dtype=int32),\n",
       " array(36, dtype=int32),\n",
       " array(204, dtype=int32),\n",
       " array(45, dtype=int32),\n",
       " array(285, dtype=int32)]"
      ]
     },
     "execution_count": 72,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "\n",
    "f_both(10)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Let's generalize this power thing a little awkwardly. Here is a list of powers"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 76,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "powers = [0, 1, 2, 3]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Use these powers as `non_sequences` and write a reduce function that gives back the accumulated sums of all these powers. Hint: use `*args` (never do that in practice) and figure out where which input arguments are placed."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 79,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "# def accumulate_all_powers(...)\n",
    "\n",
    "#f_all_powers = ..."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 80,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "# %load scan_tutorial_solutions/sol_theano_reduce_awkward_multiple_powers.py\n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 78,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "[array(10, dtype=int32),\n",
       " array(45, dtype=int32),\n",
       " array(285, dtype=int32),\n",
       " array(2025, dtype=int32)]"
      ]
     },
     "execution_count": 78,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "\n",
    "f_all_powers(10)\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Figure out a way to do this less awkwardly by vectorizing the power operation. Now the powers can even be symbolic:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 81,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "powvec = T.ivector()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Write a function that gives the same output that is much less awkward"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 83,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "# def accumulate_powers(...)\n",
    "\n",
    "# f_acc_powers = theano.function(...)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 84,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "# %load scan_tutorial_solutions/sol_theano_reduce_vectorized_multiple_powers.py\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 85,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([    10,     45,    285,   2025, 120825], dtype=int32)"
      ]
     },
     "execution_count": 85,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "f_acc_powers(10, np.array([0, 1, 2, 3, 5]).astype('int32'))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now recall the function that outputs a list of monomials from above. Reduce it to return a polynomial."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 89,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "# ...\n",
    "\n",
    "# f_poly = theano.function ..."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 90,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "# %load scan_tutorial_solutions/sol_theano_reduce_monomials_from_map.py\n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 88,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array(256.0)"
      ]
     },
     "execution_count": 88,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "\n",
    "f_poly(3.)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The same effect can be obtained by solely using reduce. Using a function which feeds back and accumulates `x` raised to a power as well as the polynomial expression, this can be done."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 92,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "# def accumulate_polynomial(...):\n",
    "\n",
    "# f_r_poly = theano.function(...)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 96,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "# %load scan_tutorial_solutions/sol_theano_reduce_monomials.py\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 95,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array(81.0, dtype=float32)"
      ]
     },
     "execution_count": 95,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "\n",
    "f_r_poly(2.)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "What about that Fibonacci sequence from the beginning? Can we calculate it using `theano.reduce`? Well, yes, we can, but it's **CHEATING**. Apart from simply exploiting the closed form formula, we can emulate looking back two steps into history. This is, however, a house-keeping trick and not necessarily the way we want to solve this problem in the end."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 97,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array(13, dtype=int32)"
      ]
     },
     "execution_count": 97,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "def fib_acc(i, old_value, older_value):\n",
    "    return old_value + older_value, old_value\n",
    "\n",
    "(fib_expr, _), updates = theano.reduce(fib_acc,\n",
    "                                        sequences=T.arange(n),\n",
    "                                        outputs_info=[np.int32(1), np.int32(0)])\n",
    "\n",
    "f_fib = theano.function([n], fib_expr)\n",
    "f_fib(6)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Can we output an accumulative sequence, as does our function `accumulate` from above? Well, 'unfortunately', yes, we can, but it's **CHEATING**, too! Take as an example accumulating the squares of a sequence:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 53,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([  0.,   1.,   3.,   6.,  10.,  15.,  21.,  28.,  36.,  45.], dtype=float32)"
      ]
     },
     "execution_count": 53,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "seq = T.fvector()\n",
    "\n",
    "def acc_squares(u, i, acc_array):\n",
    "    return T.inc_subtensor(acc_array[i:], u)\n",
    "\n",
    "acc_expr, updates = theano.reduce(acc_squares,\n",
    "                                  sequences=[seq, T.arange(seq.shape[0])],\n",
    "                                  outputs_info=T.zeros_like(seq))\n",
    "\n",
    "f_acc_expr = theano.function([seq], acc_expr)\n",
    "\n",
    "f_acc_expr(np.arange(10.).astype('float32'))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "You may ask: Why exactly is this cheating? The correct answer is that `theano.reduce` was not made for these settings and will be inefficient if misused. Especially the last example shows this: There is no need to increment a large part of a tensor at each iteration in order to achieve what is asked.\n",
    "\n",
    "Hence why we shall now move on to the most general of theano looping functions: `theano.scan`\n",
    "\n",
    "## Generalizing to generalized generality: The generalization of `theano.scan` is `theano.scan`\n",
    "\n",
    "With `theano.scan` you can do everything. Contrary to `reduce` it natively outputs all sequence values. It can also provide access to points not in the immediate history of the sequences, permitting arbitrary taps as long as the initialization is long enough, too. The general setup with respect to `sequences`, `outputs_info`, `non_sequences` is exactly the same as in `theano.reduce`.\n",
    "\n",
    "Let's start by once again adding up integers:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 54,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([ 0,  1,  3,  6, 10, 15, 21, 28, 36, 45], dtype=int32)"
      ]
     },
     "execution_count": 54,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "scan_ints_expr, updates = theano.scan(add,\n",
    "                                      sequences=T.arange(n),\n",
    "                                      outputs_info=np.int32(0))\n",
    "f_scan_integers = theano.function([n], scan_ints_expr)\n",
    "\n",
    "f_scan_integers(10)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Once again, please do the same for adding squares (e.g. using `add_square_a_to_b` from above)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "# ...\n",
    "\n",
    "# f_scan_squares = theano.function(...)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 107,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "# %load scan_tutorial_solutions/sol_theano_scan_accumulate_squares.py\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 106,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([  0,   1,   5,  14,  30,  55,  91, 140, 204, 285], dtype=int32)"
      ]
     },
     "execution_count": 106,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "\n",
    "f_scan_squares(10)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Similary to with `theano.reduce` before, create a scan that outputs both sums of integers and sums of squares at the same time. Further, have it output $\\frac{n^2 + n}{2} + 1$ at every iteration (don't feed it back, place `None` in `outputs_info`), and check that this value is always 1 higher than the sum of integers. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 111,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "# def acc_ints_and_squares(...)\n",
    "\n",
    "# f_scan_ints_and_squares = theano.function(...)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 112,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "# %load scan_tutorial_solutions/sol_theano_scan_accumulate_ints_and_squares.py\n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 110,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "[array([ 1,  2,  4,  7, 11, 16, 22, 29, 37, 46], dtype=int32),\n",
       " array([ 0,  1,  3,  6, 10, 15, 21, 28, 36, 45], dtype=int32),\n",
       " array([  0,   1,   5,  14,  30,  55,  91, 140, 204, 285], dtype=int32)]"
      ]
     },
     "execution_count": 110,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "\n",
    "f_scan_ints_and_squares(10)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Let us tackle the concept of `taps` by finally addressing the fibonacci sequence without resorting to tricks. In order to do so, take a look at the following `scan`."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 122,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([ 9.,  1.,  9.,  1.,  9.,  1.,  9.,  1.,  9.,  1.], dtype=float32)"
      ]
     },
     "execution_count": 122,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "def identity(x):\n",
    "    return x\n",
    "\n",
    "alternating_expr, updates = theano.scan(\n",
    "        identity,\n",
    "        sequences=None,\n",
    "        outputs_info=[dict(initial=np.float32([9., 1.]), taps=[-2])],\n",
    "        n_steps=10)\n",
    "\n",
    "f_alternating_expr = theano.function([], alternating_expr)\n",
    "\n",
    "f_alternating_expr()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Take a close look at the above function. We have introduced a more general `outputs_info`: It takes initial values, but also `taps`, indicating where in history to look and which values to give to the scanning function `identity`. `identity` lives up to its name and just outputs its input. We observe an oscillation due to the lag.\n",
    "\n",
    "One can give several taps to the function:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 126,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([ 0., -1., -1.,  0.,  1.,  1.,  0., -1., -1.,  0.], dtype=float32)"
      ]
     },
     "execution_count": 126,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "def minus(x, y):\n",
    "    return x - y\n",
    "\n",
    "another_alternating_expr, updates = theano.scan(\n",
    "        minus,\n",
    "        sequences=None,\n",
    "        outputs_info=[dict(initial=np.float32([1., 1.]), taps=[-1, -2])],\n",
    "        n_steps=10)\n",
    "\n",
    "f_alternating2 = theano.function([], another_alternating_expr)\n",
    "f_alternating2()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Using two taps and appropriate initial values, write the fibonacci series up to `n`."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 133,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "# def fib_acc(...)\n",
    "\n",
    "# f_fib_scan = theano.function ..."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 132,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "# %load scan_tutorial_solutions/sol_theano_scan_fibonacci.py\n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 130,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([ 1,  2,  3,  5,  8, 13, 21, 34, 55, 89], dtype=int32)"
      ]
     },
     "execution_count": 130,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "f_fib_scan(10)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Bonus: Write an awkward convolution function using this mechanism"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### An oscillating dynamical system"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We will use `theano.scan` to create a small sinusoidal oscillating dynamical system by iterative multiplication of a small rotation generator. Let us define an angle increment and its corresponding 2D rotation matrix:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 146,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "theta = 2 * np.pi / 720\n",
    "rotation = np.array([[np.cos(theta), -np.sin(theta)],\n",
    "                    [np.sin(theta), np.cos(theta)]])\n",
    "a = 1.000\n",
    "rot = T.constant(a * rotation.astype('float32'))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Using an initial vector `[1., 0.]`, apply this rotation generator iteratively for `n` iterations, keeping all outputs."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 151,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "# def rotation_increment(old_value):\n",
    "\n",
    "# f_oscillation = theano.function(...)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 152,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "# %load scan_tutorial_solutions/sol_theano_scan_oscillation.py\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 153,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "[<matplotlib.lines.Line2D at 0x7fcafd637190>,\n",
       " <matplotlib.lines.Line2D at 0x7fcafd637410>]"
      ]
     },
     "execution_count": 153,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYMAAAEACAYAAABRQBpkAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAIABJREFUeJztnXm4VMWZ/z9vQIw74IYI7rhvgAIK6NW44L4mauIes02M\nyS+ZiclkMoczM0k0mRmNMSYmE42aaEyMGncl0asIgiKiuKAgoIKK+77L+/uj6krTdN/byzmn3u5b\nn+fh4d7u06e+t9+q89by1luiqkQikUikd/Op0AIikUgkEp7oDCKRSCQSnUEkEolEojOIRCKRCNEZ\nRCKRSIToDCKRSCRCBs5ARC4SkSUiMruba84Tkbki8qCIDG+2zEgkEolkSxYjg4uBCdXeFJEDgS1U\ndRjwZeBXGZQZiUQikQxp2hmo6mTg1W4uORS4xF87HegvIus3W24kEolEsqOINYMNgWdKfl8EDCmg\n3EgkEonUSFELyFL2e8yBEYlEIoboW0AZi4GhJb8P8a8th4hEBxGJRCINoKrlHe6GbtL0P2ATYHaV\n9w4EbvI/jwGmVblO3f+6Cei9oH8E7ZuFvh71T2RjJjKXiZzLRFaucs1hTOQFNuGqIjQ1aY+JoTUU\nrRN0Y9B5oP8LWtmG6OGgL4AeU5ROJnIME1nCRA6v8v7Kvt7NYyIbW/k+W8XuRegEPcbXmyNCa6qi\nU7O4TxahpVcAU4GtROQZETlVRL4iIl/xKm8C5ovIPOBC4J+6d04sBPYE1gF+L7LCFFOmSCqDgE7g\nfE30W5ro+xV1Jfo3YG+GsL+kclSemiL1IcIg4E7gPFW+rUplGyrXAp8BzhXhiNx1pXIkcA6wjyZ6\nbUVNib6viX4LOA+4Q1IZnLeuSO2IcBTOhp9R5ZrQevIki2ii41R1sKr2U9WhqnqRql6oqheWXHO6\nqm6hqjup6sye78m7wOHAFsD3m9VYDUllZeBq4BJN9Oc96kr0YebwR+DXkspOeemK1I4IXTa8WJXz\nerpeldnAQcBvRMjNhr5+XAgcoolW3YPzia5EzwN+D/zV18tIcLYeBPwaONjXm7bG7A5k7xCOBL4u\nwp45FXM2sAT4j5o/8RJ/Ar6Ba7Sr5aSrWTpDC6iRzgzu8VPgOeA/a/2AKjOBM4CrRKjFhp31CPL1\n4q/ANzXR++v46H/h/paz6ymvhM4GP1c0naEF9IQIq8M5hwDf8PWl7RE/5xQcEVGtsAgiwiHAucCO\nqrydWXmp7AH8CdhBE325gc9fBryqiZ6RlaZIffhOwuXADqq80sDn/wC8rMo3M9WVyi+A/proCQ18\ndm1gNnCMJjo5S12R2hHhfGANVU4KraUnqj0768XsyKALVa4HpuB6TZngh+G/A77WiCPwfBM4UlIZ\nm5WuSO2I8GmcDb/aiCPwnAEcLcLumely9eFIaMzB+Pr4NeCiOF0UBhHG4aapvxVaS5GYdwae7wDH\ni7BVRvf7BjDHLwo3hCb6CnAmcK6k0irfYzvxDeAR31loCO9EnA0lg2AKVw/OAc709aMxXa5ePoZz\nVpEC8fXgHOC7qt1mVmg7WuIhpsqLwFnAfzd7L0llHeB7wL80ey/gCuBj4AsZ3CtSIyKsi3uIfzeD\n213u//98Bvc6DrfB8vKeLqyBfwbOlFTWy+Bekdr5PLAUN4XcqzC/ZrDsffoBc4ATVJnScDmpnAv0\n1URPb/QeZfcbC/wRGKaJfpjFPSPdI8LPAVHNpufspwUuA7ZUpSEbSir9gLnA8VnN9Usq5wEfaaLf\nzuJ+ke7xkWlPAF9Q5e7Qemql16wZdKHKB8BPgB82eg9JZX3gROBHmelKdAruIVD3YmGkfvyeghOA\nH2d1T9/wnwSOb+I2J+KmHrNc9D0LODmODgrjROCxVnIEWdIyzsBzCbCtCKMb/Pz/A67QRJ/LUBO4\nsMYfSCpFpPfo7Xwb+KMqz2d8X2dDqT9Fi7f798iwkwGgiT6Lm4r8Tpb3jayIt3vmNmwlWsoZ+NHB\n2cAP6v2spDIA+BIuLj1bXYnehcvGemzW944sQ4SBwBeBn2V9b1XuBJ4FPtfAxz8HPOvrQdacDXzJ\n199IfhwLLFKl14bztpQz8FwM7CbCZnV+7jTgRk30qRw0gVvc/qakkmv6jF7Ol4AbVHk6p/v/N3WG\nhHp7f4fGN4p1iyb6NHATcGoe94+AT3nzbXKyYavQcs5AlXdwDuHrtX5GUumDy4n0i7x04Rrs2tDw\nFFakG0Tog4u/z9OGNwLriTCqjs+MBvoDN+cjCXB/89d9PY5kzxhgTeCW0EJC0nLOwHMBcFKNqQTA\nZU59QRO9Ly9BmujHwC9x8e+R7DkIeF6VGXkVoEojNjwd+KUmujQfVaCJTgdewn0HkexxNlRys2Er\n0JLOwGc2vZvaoz++Dpyfm6BlXAQcKKlsUEBZvY3TKc6GB4vQ49GsPjrtINxINW9+QexoZI6384G4\nJIG9mpZ0Bp4LgC/3dJGksjkwAvhL3oI00VdxCcpOzLus3oQIWwA7UYQN3a7kq6GmnDSnAn/1ds+b\nPwM7Syr1rpVFuueLwFW9bbdxJVrZGfwDN7+7Yw/XnQRcrom+V4AmcL3EU+JCcqY4G1Y5pyAHLgZO\n7u4sDW/fk3H5kXLHn7NxBbU5qUgNePueTEE2tE7LOgM/v3sJzpgV8bliTvTXFcVU3PcaF5IzwOeK\nKdqGU4CVoNuF5NG41BPTClHkuBg4KebCyowxuPPYp4cWYoFWr1SXAF8QYaUq748H3gBmFSVIE1Xc\n/OMpRZXZ5uwBvKZaoA2VLhue3M1lJ+IORSoyn8ss4HXI7XyP3oazobN3r6elnYEqc3G5RA6scslJ\nwKUFN1iAS4HPSiqrFFxuO3ISxY4KurgU+JwIK9jQp5b+HPCHIgX5enwxsaPRND4FeuE2tExLOwPP\nZVTIOOlPmzoCl0SuUDTRRcBMqjupSA340OHDCWFD5RlcT/yACm8fAjyU4wbG7rgcODR2NJrmEGBW\njhsYW452cAZXAxMq7Dk4GLg3hzxEtXIlcEygstuFQ4BpqiwJVH41G36eQD1KTfQFYAaxo9Esn8d1\nJCOelncGqryEW8Qr35BzNC4cLxRXA/sbPie5FbBgw+U6GpLK6sA+wLXBVDkn1UgOpQggwhrA3kDD\nh1u1Iy3vDDzL9eD8A3g/AjZYf3zhPbgRSqRO/AN4XwI22CodjYOAqc2cZJYB1wATYkejYQ4CpsS9\nBcvTLs7gWmAf7/HBzfNOb+J846yIPbjGORC4p4nzjbPizyxvw6MpYPNbd2iiL+E6GjE9RWMEt6FF\n2sIZ+AfG3bg5ZnDGviqcok9wTiqVNXq8MlKOFRteA+wrwholI04L0wtxTaoBLIw4rdIWzsDzF+Bo\nH2UxgbBzusAn6SnuJvbg6kKEVYH9sWBD19GYihupHIALSngprCrAfTf7SiqrhhbSYhyIC0oIPeI0\nRzs5gxuBz/Bu/0OA+33UhQX+xrIRS6Q2JgD3+Tl7C3TZ0Mz0gu9o3IdbzI7UjhkbWqNtnIEqLwIP\n8c46X8ZFgVjhBuAASaXaLunIihyBLRteT5/3DkA5AFvTC9cBh4YW0Sr4A+8nYMuGZmgbZwBAn/ev\nY/Xnx+IewCbw59jOw6XGiPSAP8TmACzZUFnMsJtf4r0BizTRUHseKnE9cEg89KZm9gQe9h3HSBnt\n5QyOOOFp3hjal4lqbVdh7MHVzhjcWbTPhBayHLv+6hUeP9RUKKImOh94ge4T6kWWcTCGOhnWaC9n\nsP1fduLJfV/HnV9gietwKQRiWuueORi3/mMGSUXY+K4h3Pv1jbpLax2I2NGoAW+36Ay6ob2cARzM\n4tG3AIeFFlLGbFy64+1DC2kBLDbYrejzPjw7cimwbWgxZfwNe/XdIlsDfYCHQwuxSts4A0llI2AD\nHj/kQoz1lHy2ydiD6wERNgYGAfeG1lLGwQg3wKcs2nAGMEBS2SK0EOMcDNwQ01VXp22cAS6W/xY+\nWOMeYKgIg0MLKuMmXCRDpDoHATf5g4ss0TV1Zc6GmuhS4GaM6TKIuelHa7STM+jy/B/hjsTcL7Ce\ncu7CnWG7VmghhjHXYCWVAcBI4HacDUeIsGZYVStwK9EZVEWEgcBw4I7QWizTFs7A78LcA9co8P/v\nH07Rimii7+J2sn4mtBaL+F3H44HbQmspY3/gLk30HVXeweUE2juwpnL+DuzhD92JrMj+wJ2qvBta\niGXawhngHiIPaqKv+d9vxeWTsRZ/fQvGnJQh9sAdNvJaj1cWy344u3VxC8Z64T4h46PA2NBajFJu\nw0gF2sUZ7EtJj1KVRcDzuOG9JW7FnXFgLTzRAsvZ0ALeTuW6nA3thZiaGw1bwNtpX2BSaC3WaSdn\nUG5si43jMVx421ahhRjEYoPdClDcOdtdPAqsBAwLoqg6cd2gMlsDHwNzQwuxTss7A0llELARLmlX\nKRaH80qcKloBETYAhuLCJC2xHzDJ2w0AH5porm7hwnGHSiobhBZiDGfDGFLaIy3vDHBZGzs10Y/K\nXp8M7CBC/wCausPiiCU0+wB3+EgwS1QbrZizoa//FqPoQmNu+tEq7eAMKjZYVd4DpmAveucfwHhJ\n5dOhhRjC3BSRzzK7B85e5fwdGO+zYFoiThWVIEI/XHBJJRtGymhpZ1Blga+USRhzBj4P/SO4hGy9\nHr/Atw/GnAHOPk9qoitkuPRn5z6GPRtOAvaOAQqfMAaYq0ro429bgpZ2Brg8Me8DT1Z5/3bsxYSD\nXV0h2A54j+o2DEVP0wvmbKiJPgW8hb38SaEwN+K0TKs7g30pW+Ar4yFgHYOpKe7A2IMkIM6G9hb4\nenqQWLWhVV0hiM6gDtrCGVR7U5WlwJ3AXoUpqo0puNQUq4UWYgBzDVZS6Y/LMDulm8umAMP9AeuW\nuB179b1wfODIdrhd/5EaaFlnIKn0BcbRc74Ri8P5d4D7cfp7LSJ02bAzsJRyxgPTNdH3ql2gytvA\nTOzt+r0D6Iinn7EH7uD7qjaMLE/LOgNc4qmnNdGeDk035ww8VnUVyQhgoaGD77vYk9oclDkbaqLP\n4Xbf7xRaS2BqtWHE07QzEJEJIjJHROaKyJkV3u8QkddF5AH/79+aLdPTQW3GngN8WoRNMyo3K8w9\nSALQgc0G20GLOgOPVV1F0oHNumWWppyBiPQBzsfFNm8LHCci21S49E5VHe7//VczZZbQQQ3G9guT\nFhvHdGArPz/dW+nAremYwacY35oVd7RXYjqwjQjW0pJbrO+F4dcLtqQ2G0Y8zY4MRgHzVHWhqn4I\n/InKR/BlGvfs1wvG4vLL18IdGFtU00Q/AKbh5jZ7HX69oB4bFsU43HrB+z1dqMr72LThncBYv3Gu\nN+JsqHwQWkgr0awz2BB4puT3Rf61UhTYXUQeFJGbRCSLGOidgUWVNgRV4XZgb4OZJntzD86t+Si1\n2rAoOqhvesGcDX1K6/nALqG1BKKDOEVUN32b/HwtseEzgaGq+o6IHABcixvCrYCITCz5tVNVO6vc\ns4P6jL0A+BCXhXJOHZ/LmzuA34YWEYgObDbYDuDbdVx/B3BhPlKaomu/wT2hhQSgA/hWaBF5ISId\nuL8xU5p1Botx2Sa7GIobHXyCqr5Z8vPNInKBiAxU1VfKb6aqE2sstwP4fa0iVVER7sKFDFpyBjOB\nTSWVAT5NRW+iA7g4tIhS/HrBNrgMoLVyP7CZCP2NHcxzF/DV0CKKxq/f1Lrm05L4TnJn1+8ikmRx\n32aniWYAw0RkExHpBxwDXFd6gYisL+JypYjIKEAqOYJa8fHT46h/rrnLGZhBE/0QN+dsLVY9Vwyv\nF4wF7q1lvaALVT7ELSRbs+HdwG69cL9B13pBzTaMOJpyBqr6EXA6Llvio8CVqvqYiHxFRL7iLzsa\nmC0is4BzgWObKRO3XrBYE32hzs9Nxpgz8FjVlSduzUep14Z500FjU1fmbOj33yyi9+036MDm9KN5\nmp0mQlVvBm4ue+3Ckp9/Cfyy2XJK6KAxYz8OrCbCUNXlFr1DMxn4cWgRBdOBsZBSTwfwLw18bjKQ\nVch0lkzGRTrNDC2kQDqA74QW0Yq04g7kPWngQeL3G5jrweGmGHaQVFYNLaRAzO0OlVTWxO2Vmd7A\nx6cBO4pgzYbmpkbzRIQ1qX/NJ+JpKWcgqXwKV7kbnWvu6imZwecpegh7ufFzQYRG13zyZndgRnf5\niKqhyjvAbGB05qqaYzLuICVrIdV5MRa4L+YjaoyWcga4LIQvaaLPN/h5qz0lq7ryYDvgRVWWhBZS\nxjjcomujmLOhJvoM8DYupLo30KwNezWt5gzG0pyxHwSGiLBORnqywuL0VV40a8O8aFaXVRta1ZUH\nVutWS9CKzqC7HPPdosrHuE041lJHTwFG95L0AeNowoZ54L/3XWhug9YUYIwI1mxobmo0D/x5x7vg\n1m8iDdBqziCLYaDF4fxruGMfR4TWUgBNOfScGA4s0ERfb/QG/lzk+f5eljBX33NiOPCkKg3bsLfT\nMs5AUhkMrIELEW0Gq8Nmq7oyQ4QhwOo0b8OsyWp6waINnwBWkVQ2Ci0kZ+IUUZO0jDPAGXtqN+cd\n18p9wLYirJ6BpizpDcP5scAUg+cdZzVaMWdD317uxp6Tyhpz04+tRis5g0wiBXzY2QPAbk0rypbJ\nwDgfPtuumOu9+bDLrKJQnA3tZce1OGLJDP99m6tbrUYrPXiynGuegjFn4I8rfI0qGV3bBIu9t82A\nj4Cnm72RKs8Cb2DPhubqe8ZsAXyg2rwNezMt4QwkldVxmQhnZHTLqbhNRtawqqtpRFgDF+9+f2gt\nZYwDpmQw/diFRRvOAjb3u6zbEYtBCS1HSzgD3M7OWfVkk+yBe3BhgNb+fosPkqwYDcw0mE0y6+kF\nczb02XHvx94O6ayIU0QZYO1hWI1MPb8/XetFXC4aS9yDsQdJhlicIoLse5VWbWhVVxZYrVstRas4\ngzy2mZvrweHy2wyVVAaGFpID5npv/nseissNlRUPARuJMCDDe2aBxfreND6bwGBc24k0gXln4A/n\nGI2rzFkyFWOLaproR7iMi201nPeH2eRhw2bZHZjuv/dMUOUjXPiyNRveg9vlbr7N14mzofveI03Q\nChVjB9xhNi9nfF+rPSWrupphR+AZVRo+4S4n8hqtmLOhJvoi8AL2pkabxdyIs1VpBWeQ13zgo8D6\nBpPWmXuQZIDVOd28dFm1oVVdzWC1brUcreAMcgkb80nr7sXYVBEu0dYoSaXpU+gMYS70T1Lph8tn\n08hhNj3hbCjNnySYMW21iCzCyrgjVPOwYa+jFZzBGJrLJtkd5npKmuirwDO46bF2IU8bNsrOwDxN\n9M2sb+ynwxYD22d97yYxV9+bZDjwhCpvhRbSDph2BpLKIGAtYG5ORVhtHFZ11Y0Ig3HJ6fKyYaOM\nId90xxZt6KZGU1k3tJCMyNuGvQrTzgAXkTFdE12a0/2nAyMM5qA3F+nUBM6G9pLTFeEMTNlQE/0Y\n9ze3yxGr0RlkSCs4g9yM7XOfLwB2yquMBrHYq2wUqw22N44MwK6uRrBat1oS686gCGNbXFR7Augv\nqWwQWkgGmGuwksr6wADc95wXjwMDRBiUYxmNYLG+143/XtfE3vRjy2LWGfjNZrvgIn7yxFxPyU+L\n3YOxaYZ68dE0I8jfhvWS9/QjqizFOUFrNpwOjGyDI1a7ph9zs2Fvw6wzALbDbTZ7NedyzDkDj1Vd\n9bA9brPZa6GFlFHUaMWcDf3RnvNx0VStjLkRZ6tj2RkUZex5wCr+SEZLmHuQNIDVBttrnYHHqq56\nsFq3WhbrziD3zSQ+yuUe7EVY3AfsJKmsHFpIE5hrsAVOP+LL2FmEfgWUVQ8W63vN+OnHomzYa7Du\nDIp6kEzHWOPQRN/CjVqsRTrVgzlngJt+fFYTzT1Pkt8M9ST2bGiuvtfJdsAiVfKeQu5VmHQGkkp/\nXGrhhwsqchr2skyCXV09IsJAXGrhR0JrKaNoB2XRhl3RauuFFtIgFjsZLY9JZwCMAu7PMrVwD9wH\nDDe4+Ww69h4ktTIKmOFzQFmi6AeJORv6KKpWTpUenUEOWHUGhRpblTeAp7CXD6iVh/NWG2wIZ2DR\nhlZ11YLVutXSRGewDIvD+ceAdSUVa2m2a8Fcgw0w/QjOhuuJsHaBZdaCxfreI/4EuSHYm35secw5\nA0lF8BtKCi7aXE/JD+ctnprVLSJ8ijA27Imipx+7UqXPwJ4N7wV29dFVrUTX9GM82SxjzDkDYAvg\nTU30uYLLtdpTsqqrO7YEXlVlSWghZYQarZizoSb6Eu7ks61Da6kTcyPOdsGiMyhkf0EFHgEGGzzI\n3NwCZA1YbbChdFm1oVVd3WG1brU8Vp1B4cb2w/mZuGGoJabjTj6zaKtqmGuwAacf8WWO9tNnljA3\nNdodhqcf2wJrlRPCPkgsDudfAF7FTb20CuacAeGmH/HTZa8Bw4ouuwfM1fceGAa8rsrzoYW0I6ac\ngaSyKm4O84FAEqz2lKzqWgERVsc12lmhtZQR2kFZtOGDwBaSyuqhhdRIaBu2NaacATASeFgTfS9Q\n+V3DeQlUfjVaqQe3C/CgKh+EFlJG6AeJORtqoh/gHMIuobXUSGgbtjXWnEFQY6vyLPA2bkrBEhZ7\nldWw2mBD67JqQ6u6KhHahm1NdAYrYq4Hh5s229JPo1nHgg2Xw8D0I77srUSwZkOL9X0FRFgNt25m\nbfqxbbDmDHYj/IPEXLidJvo+btfsyNBausNPr5lzBrjvbbb/HoOgynu48OURoTRUYTowxkdbWWYX\n4CFVgtmw3bHmDPoCCwNrsDpstqqrlI2BpcAzoYWUYcVBWbThU7jnwNDQQnrAig3bFmvOYJomqoE1\nzAS2FWGVwDrKaYXh/Bhgmj8wyBJWHiTmbOjbmzldFbBiw7bFnDMILUCVd4A5wPDQWsqw2KssZzQG\nbFiKn/6wMP0Idm1oVRdgevqxrYjOoDIWe0rzgZUllQ1DC+kGiw12CK6ePxVaCO7Us1VEGBxaSBkW\n63spXVNYTwdV0eZYcwb3hRbgsbiIrBjU1YUIKwM7AveH1lKGm7oKP/3Ydd62RRvOAHaWVKwd7tSF\n1enHtqJpZyAiE0RkjojMFZEzq1xznn//QRGpOv2iib7ZrJ6MmIbNYbNVXQA7A0+o8nZoIWWESnxY\nDXM21ETfABbgnLlFrNmwLWnKGYhIH+B8YAKwLXCciGxTds2BwBaqOgz4MvCrZsosiLnAmiIMCi2k\nDIu9yi4sThGBPV1WbWhVF9izYVvS7MhgFDBPVReq6ofAn4DDyq45FLgEQFWnA/1FZP0my80VPxy1\neEbsfcBISaVvaCEVMJdNUlLpB+yEmwaxwr3ASBGs2dDkIrIIFm3YljTrDDZk+ZjyRf61nq4Z0mS5\nRWBuUU0TfQ23iLZ9aC0VsDiU3xFYYGj6EVVew7WB7UJrKcNcfffsBDypyluhhZQiqewqqVibOWiK\nZnsntS7olO9urPg5EZlY8munqnY2oCkrpgP/HLD8anQN581syxdhPWAg8HhoLWWYC3X1dNnwwdBC\nSngU2EBSGaiJvhJaTAlWbfgz4MdQfDptEekAOrK+b7POYDHL71wciuv1dHfNEP/aCqjqxCb1ZIk7\nI1bo4w++sULXAuSFoYWUMBq4V5WloYWUMQa4M7SICnTZ8DehhXShiX4sqczATf3eElpPCWOAO0KL\nKMVP047EPSMKx3eSOz/RI5Jkcd9mp4lmAMNEZBMR6QccA1xXds11wIkAIjIGeE1VrZ2NuwKqvIzz\n+tv0dG3BWFzos9p7s6rLog3Bpi6LNtwOWOSnbduGppyBqn4EnA7cihtmXqmqj4nIV0TkK/6am4D5\nIjIP15v9pyY1F4nFRbWHgY0klf6hhZRgcfF4bWAQ8FhoLRWYDWwswlqhhZRhqr6LsDawHi4jgCXa\nMrqp6X0Gqnqzqm6lqluo6k/8axeq6oUl15zu399JVWc2W2aBmFtU00Q/wm3s2jW0FvjkXNpRBBoy\nd8MoYIYmammKDwBVPsLlwDJhwxKm4c7btpLBdDQww9g0LRjs/GSBtR3I1jC3QchjSdfWwIuqvBha\nSBkWpxdKsWRDADTR54E3sXNWs1UbxpFBL+QhYFMR1gwtpAxLDxKLIaVgV1cXlmxYiiVd5mzop2eH\n4qZr24roDLpBlQ9xJ1RZG85bOpDEXO9NUumaujL1ICljGjDG6HnbwZ1ByfSjNRvuCjzgp2vbiugM\nesZE4yhFE30WeAsbZzWb673hpjne8NMeJvHnbb8DbB5aSxlW6vuWwGuqWIs8NNf5yYroDHrGSuMo\nJ7guEVbHOSRLm6egdRb4gtuwAu6s5lRWC6zDqg0tdn4yITqDnonD+eqMxOa5tK3Se7Ngw+Xw50TP\nJvx52+Zs6KdlzenKiugMekCVxcB7wGahtZRh4UFitZdkVVc5FmxYCQu6LNpwM+B9TbRiBoVWJzqD\n2rDQOMp5ANhaUlk1oAZzvST/fWyN+36s8wCwjQghbViJoPXdfx9bYc+G5up7lkRnUBvmnIEm+h4u\nvC3IcL7kXFprvbcRwKOa6LuhhfSEKu8Cj+A0W2IasFvAaLWRwCOqvBeo/GpYrO+ZEZ1BbZhzBp6Q\nuobgEh0uDFR+NawuPFbDYt16CpdpeGhPF+aEVRvGkUGEmcC2IqwSWkgZIR8ko7F5Lm2rNVhzzsCf\nFx28bgUquyKSyqdx54hYO+M7M6IzqAE/nH+UOJwvxeqQ2aquaphzBp6QuizacGfgcU30ndBC8iI6\ng9qx2GgXAn0Ic3Kcxd7bYGA1YF5oLXWwAOgnYu70vyD1XYQNgVWAJ4suuwcsOqhMic6gdsw5g1DD\neRFWAobjzmS2xGhguv9eWgI/zWaubuHOKtlJUlm54HKdDeP0Y+FEZ1A7FhsshNG1A7BQlTcKLrcn\nWrXBmqtbmuhbwFzcGcRFYtWGcWQQ+YT5wKfjcB6w2zCs6uoJc87AE+sWIKmsD/QHngitJU+iM6iR\nkuG8qcNucMP5nSWVfgWWaS70T1LpQ8BzaZvkPmBnEYq0YS0U6gxE6IsL0rBmQ3fGd6LWzvjOlOgM\n6sNcD04BXKggAAAbK0lEQVQTfRO3YFrkcN7iUH47YLEm+mpoIfWiypu4keeOobWUUXR9d2cLK9bO\nFrZY3zMnOoP6MOcMPIXpEmEAsCFu56wlzE0v1InFuvUEMMBPkxSBVRta1ZUp0RnUx33AcB9NY4ki\nHyTubGGb59K2cu/NnDPw0yLTKW5q1JwN/fTjLkRnECnFR88soHcP5632kqzqqhVzzsDT2+vWNsAL\nmujLoYXkTXQG9WOx0T4OrC2prFdAWbsDUwsop2YklQHARrgzq1uVOcA6IqwbWkgZhdR3EQbiciHN\nzrusOjFX3/MiOoP6MecMihrOi9AH97ffk2c5DTAGuK+Vz6VVZSkuisZatNq9wC5+uiRPxgD3qmLN\nhtEZRKpizhl4itC1LbBElRdzLqde2qXBmqtbmugrwLO4SJ88sWpDq7oyJzqD+nkMWE+EdUILKaOI\nB4nVhmFVV72YcwaeXlm3JJV1gfVxSSrbnugM6sT4cH7XnIfzu2GvwfYFdsVYFEqDTAd29dNxlsjV\nGfjNZrtgz4ZjcLmurEXO5UJ0Bo1hrgfnox2ew03l5IW53hsux/xiP53R0qjyMrAEF8Fiibzr+w7A\nM6pY2zBosb7nRnQGjWHOGXhy0+WjXNbD3pC53Rqsxbr1MDBUUumf0/2t2tCqrlyIzqAxpgOjjA7n\nd8vp3rvhUgtby8/Sbg02Txs2hI/Sup/8pkZ3x1iEmqSyEi7XlbV9D7kRnUEDqPIS8Dz5R1jUy1Rc\nw8oDqw9dq7oaJU8bNkNvq1s7A/M10ddDCymK6AwaZwowLrSIMh4GNvBREFljrsFKKhvgUgs/HlpL\nhswGBhuMVsulvoswGFgTe+mhzdX3vInOoHHuBsaGFlGKj3q4h4x7cD618gjsDZl3A6a1U2phn/Np\nOvZGB1OBUX76JEt2A+4xOv1oauoqb6IzaJwpGHMGnjx6cDsBTxo82cxcqGtGWOxovAo8Rfap0q3a\n0Kqu3IjOoHGeAFYzePJZHk7K6pDZqq5msTgFCb2kbkkqQ4FVcOeE9BqiM2gQf/KZxdHBdNxB5qtk\neE+LDXZl3CKftVOxsmA67uSzT4cWUkamIxb/9+2ESw1vid2AqZqohhZSJNEZNIc5Z6CJvo3bC7BL\nhrc15wxwaxiP+4Pb2wpV3sKlPRkZWksZU4BxkopkdL8RwGOqvJ3R/bLCYn3PnegMmuNubA7nM+vB\niTAUWBl3LKMl2r3BWpwqWuD/3ySj+1ldpG33ulWR6AyaYyawpQhrhBZSRpYPkt2BqX5azBJWHyRZ\nYXHU2TU1mmndyuhemSCprIpLcTIjtJaiic6gCVR5H3gAe0nrpgC7SypZ2Hecv58Z/DSFOV0Z42wo\n5tpoJk5KBKs2HAXM1kTfDS2kaKxVtFbE3FSRJvoc8CqwdQa32wO4K4P7ZMmWwPua6MLQQvJClcXA\nm8BWobWUkdUU5FbAO6o8ncG9ssRifS+E6Ayax9xw3tP0cF6E/sBmuOkwS4wHJocWUQAW69aDwCb+\nqNFmsGpDq7pyJzqD5pkKjPY52S2RRQ9uLO4owg8z0JMl4+kdvTeLm88+xIXzNptMz5wN/e7q0bjv\nvdcRnUGTqPIK8AywY2gtZWSx0Gd1yLwHvaP3ZjGiCLKrW9ZsOBx4qh3OxmiE6AyyweJw/jFggKQy\nqIl7mBsySypDgNVxf1+78wiwrgjrhxZSRlP13Ycrr4q9BIPmRitFEp1BNlhcRF6Km8JqSJcIq+J2\nh1o7inA8cHdv2B3qk7dNxV5H4x5gpKTSr8HPjwcmGwxXNtf5KZLoDLJhMjDeh8tZYjKugjfCaOAh\nVd7JUE8WWJxeyJNmbJgLmugbuNxcje5yN2dDH4YdnUGkaRYCHwLDAusopxPoaPCzVhtGbxvKd9K4\nDfOkk+bqljUbbgO8pokuDi0kFA07AxEZKCKTROQJEblNpPL5qCKyUEQeEpEHRKQdk4p1Ja27E3uN\ndiawqaSydgOfNecM/N+xETArtJYCmQFsLkKzoZxZ01B994f2DMGFqFrCXH0vmmZGBt8DJqnqlsA/\n/O+VUKBDVYer6qgmyrNOJ7BnaBGl+DDAqdQ5zSDCSrhD2a3tDh2HO8zmo9BCisKH9U7D2FQR7sE5\npoHDbsbhDrP5OAdNzWBu6qpomnEGhwKX+J8vAQ7v5lprc+l50Al0GFw36KT+HtxwYIEPm7WExemF\nIujE2KjTh18+Sf3rBuZs6NObWA2jLoxmnMH6qrrE/7wEqoa/KfB3EZkhIl9qojzrLAA+BrYILaSM\nTup/kFjtJVnVlTedGHMGnk7ao25tAvSllx1mU063u2ZFZBJQKU79B6W/qKqKSLUwsbGq+pyIrAtM\nEpE5qlqxMojIxJJfO1W1szt9llBFRT6ZR50bWE4p9wObSSoD69hMsydwWY6a6kZSWRPYlvY8zKYn\nZgDDRBigyquhxZRwJ/A14Ce1XCzCWrh8WdYOs9kTuLNVwpVFpIMcOgfdOgNV3bfaeyKyREQGqerz\nIrIB8EKVezzn/39RRK7BZQWs6AxUdWKtwo3SCewF/Dawjk/QRD+UVO7BDc//1tP1Pq3GHsAX89ZW\nJ+OBe3tjNklVPhDhHtx8+/Wh9ZRwF3CppLKSX5/qifHAdFXey1lXvewN3B5aRK34TnJn1+8ikmRx\n32amia4DTvI/nwRcW36BiKwqImv4n1cD9gNmN1GmdTpp/XWDkcBTqpWde0BaqsHmgLloNT/SnE/t\nJ7KZs6FfLzCnKwTNOIOzgH1F5Ancl3kWgIgMFpEb/TWDgMkiMgt3rusNqnpbM4KNMx+3RrJ5aCFl\ndFL7g8Rqw7Cqqyg6MeYMPJ20dt0ahmuzvXq9AJpwBqr6iqruo6pbqup+qvqaf/1ZVT3I/zxfVXf2\n/7ZX1ZrmFluVkv0GpkJM6YpVry3tsLkG6/cXbI69ueYiuQ93ql7F/TwBqam++/0Fm2LvBLG9gdtb\nZb0gT+IO5OzpxFgPzs/n9hirLsLKuP0F1kLsOoApNc5LtyWqfICzoakcWLi6MraG/QYdwN2qWNsj\nsjdwR2gRFojOIHs6gb2Mrhvs1cM1o4E5qryWv5y62Atjo5VAdNKzDQtFE30ZF1bd07qBORv6fEQd\nRGcARGeQB0/i9htYO67wH8A+PVxjborIY1VX0dRiwxC0at3aDnhDE30qtBALRGeQMX7dYBIucsoS\nM4ANJZUNurnGXIP1egfRu/IRVWMGsJFIxb0/IZkEdBOGzmBgPezlIzJX30MSnUE+3EY3jSMEmujH\nuIpfUZcIqwEjsHfk3164DUHWctkUjp9vvx17o4O7gBGSujDyCuwFdPrzGSwRnUEJ0Rnkwz+APXzC\nN0t014MbB8xS5e0C9dTCZ4hzuqV02wsPgSb6Ni7aqVpUkTkbSipdmys7A0sxQ3QGOaDKy7iUFGNC\nayljErCP32hTzn64EY0ZvE5zugIzCdjXYIBCRSflde6PPRuOBp7WRJ8PLcQK0Rnkh8Ue3HzgHWD7\nCm9PAG4pVlGPbItbjLd2Vm4wVHkSeA/33ViiWn3fDngfW/m6wDkoa/U9KNEZ5Ic5Z+BZQZc/oHx9\nXFI7S+wP3Bo3BK2Axbr1ALCepDKk7HVnQ3vnHe8P3BpahCWiM8iPKcD2Bk+oqvQg2Q+YZPDAkdh7\nq4y5aLWSAIXyxe0JGHvoSirr4I65tHZ4U1CiM8gJVd7HVTZTm4RwDXaspLJyyWsWG+yqwO7EaI9K\n3A6M8zvGLbFcR0OEVXHrZtZsuA8uQu390EIsEZ1BvtyGvR7cq8CjwFj4JGX1Z7C3wLcHMEsTfT20\nEGv4E+gewzlLS9yGC1Doeq7sCcxU5Y2AmioRp4gqEJ1BvtwEHGgw8uMm4ED/867AM6o8G1BPJeIU\nUffczDIbmsDv5H2JZakpLI44u6KbYt0qIzqDfHkc+BDYIbSQMm4EDvY/m2uwHqu6rHADy2xoidK6\nZbEHvgPwriba61NWlxOdQY74CAqLjfYBYE1JZRgGH7qSysbAOsDM0FoMMxPoL2LuzO0bgINF2AQY\niKtrljBX360QnUH+mHMGmuhS4Ebe2PBYYEvsHVB+CHCj1xmpgE/tcCNwUGgtZUwFNmXDaccDNxpM\nQXEIto4ONUN0BvlzF7CdCOuGFlLGDSztexwuBvyD0GLKOBR3rGqke0qnZEygiX4E3MKmd3weYzaU\nVNYFdsRYagwrRGeQMz7E9O/AAaG1lPEPVnthGAPnTgotpBRJZS1cOKK16CaL/B0YI8KaoYUsxxuD\nb2ejKcOwZ8ODgEma6HuhhVgkOoNiuAFrw/mJqjw1TjhxX2unh00AJmuib4UWYh1V3sRNy9jajfzH\nm5VNboeJYm2KKI44uyE6g2K4CdjPWBbTfXhm3JP0f8raprjYYOvD3JoUS3b8DO+uvQBDx79KKp/G\n7ae5KbQWq0RnUACqLMGFmVp68B7Kgr0vBw6u4fzaQvA6DsA94CK1cT0ueqdvaCEAvsMzAdErgcND\n6ylhb+BBTfSl0EKsEp1BcVwFHB1aBIAInwIO5ulxl+HOr62Wh75oxgFPaqKLQwtpFVRZCCzEjg3H\nA/NYc/FFwBH+3AALxBFnD0RnUBxXAYcb6cGNA5aoMh/4C0acFHAU8LfQIloQMx0NnA2v1UQXAE/j\nnENQJJU+uFFKrFvdEJ1BQfge3FO4nDuh+Rxwpf/5r7geXJ+Aeroa7NEs0xWpnauAI0UIa0NX/lHA\nn/1LVwGfDafoE/YAntVErZ2pYIroDIoleOPwDfZofIP1B94sIryTGk9ssA3hD7xZTPhe+B7AYlW6\nUj38FTgydEcDOIbYyeiR6AyKxTWOsD0410vS5U6esjDNcAzLepSR+jFnQ+/Yn8dnyA2BX7M4Ejcd\nGumG6AwKxPeYniVsD650iqiLq4CjQvXgfIMtnV6I1I+zYaCOhl8LO5IVbfgXwo6G9wIW+hFwpBui\nMyieK4HjQhTsG+xRlPWSfA/uWcLFhXcAT8UG2ziqPIHrhXcEkuAeusqCstf/DHw2YPjy54idjJqI\nzqB4/gB8VoRVApS9D67BVnroXgqcVLCeLj5PnNPNgkuBEwOVXdGGvqOxgACHPPmNZkcQnUFNRGdQ\nMKoswh08f2iA4k8Gfl/lvcuBQyWVNQpTA0gqq+Ma7B+KLLdNuRw4TITViyzUl3c48Mcql1xCmI7G\nobiNZk8HKLvliM4gDIU3DhEG4PL+XFHpfU30BVyG1aOK1OXLu1sTfb7gctsOv9N9MsXb8GhgsirV\nbHglsJ+kMqBATeA6PxcXXGbLEp1BGK4BdhNhUIFldqWrfrWbay6h+GmGU6g+WonUT4heeLc29Odu\n34aLNioESWVDYDfg6qLKbHWiMwiAKm8D1wLHF1jsyfTcS7oB2FFS2SRvMQCSymbAdsTDRrLkemBH\nf9JY7oiwObAtPeeTugRXB4viBOAqTfSdAstsaaIzCMdvgK/6PEG5IsLOwAZAt2cXaKLvA5cBX8lb\nk+c04HJN1NrhOi2LPz/jD8CXCyryNOCPNRyQdCswWFIZnrcgSeVTwBeJU0R1EZ1BOKYBr+MODc+b\nbwC/UuXjGq69APiij8TIDX//04Bf5llOL8XZUMjXhi4i7ovUYEN/Atqvga/nqckzAXgTuKeAstqG\n6AwCoYriGlGujUOEtXGbgX5bky4XCjiT/Od3jwVmaqJP5FxOr8PvOZhF/pu9jgVmlO1m747/w21u\nHJijJnCdn19ooppzOW1FdAZhuQIYLcJmOZZxGvA3VV6s4zPnk6OTklQE32DzKiPC+cDped1cBAHO\noA4b+oi1G3ALzvnoSmVLYCTwp7zKaFeiMwiIKu8CvwO+lcf9/UEj/0T9D92bgYGSSl5pM8YCa/ly\nIvlwE7COSG55gcYDq+HWAurhF8AZkkq/7CUBrpPxO0303Zzu37ZEZxCec4HjRVg/h3sfD8xV5f56\nPqSJfgycDfwgB034+/5ME7V2Rm7b4NeH8rehUpcNNdF7gbnAF7IWJKkM8vc9L+t79waiMwiM36hz\nOfDtLO/r8xD9APjPBm9xKbCdpLJLdqpAUhkFbE/cW1AEl+DCTEdmeVMRRgHb+Ps3wo+A7+eQGPGf\ngT9oos9lfN9eQXQGNvgpcJoI62Z4z2NxqarvbOTDPsz0Z8APM9QE8G/A2f7+kRzxYaZ52PCHwFk1\nhJNWoxN4kQyDFCSVdYFTcW0p0gDRGRhAladxeV2SLO7nQwr/I4P7/QbYSVLJ5HxdSWUcsDMuqiRS\nDL8Bhotkc3iRCOOBnYCLGr2Hj/L5AfCjDEOYf4jbs7Ioo/v1OqIzsEMKHCPC1hnc65vAbFXuaOYm\nmuh7wPeA//EbeRrGf/4c4Pv+vpEC8EEKzoZNbnD0nz8HOFOVpmyoiXbiwl/PaOY+AJLKNrh0KxOb\nvVdvJjoDI6jyMvAT4BwfttcQIgwG/gU3f5oFVwIf0ny+m5OApVRJlBfJlT8BH9N83qlTcHUhq7DN\n7wLflVQ2aPQGPkz5HOAsTfSljHT1SqIzsMX5uLQRDTVa70R+A/yyjo1A3eKH9F8DzvbJv+rX5T53\nNvDVGEFUPH6D49eAn/rOQt2IMATXWfmqv1/zutwGxwuAC/1DvRFOBtYnRhA1TXQGhvALcicD/y3C\nRg3c4mRgCC5aIzM00Vm4RvvbeqeL/PW/BS7QRB/IUlekdlR5APgV8Nt6p4v89f+H62Q8mLG0/wI2\npoEkdpLKxrgF45M10Q8z1tXraNgZiMhnReQREflYREZ0c90EEZkjInNF5MxGy+stqDILOAu4RoRV\na/2cCCNwDeOEJqI8uuPHwJrUvyj9b7gNZj/OXFGkXn4E9Af+vc7P/TuwOjnY0Ccp/ALw03rCmCWV\nVXDpqc/WRLN2UL2SZkYGs3EnVN1V7QIR6YOb+piAS3N7nIhs00SZwRGRjgKK+V/gUeAPfhdxt4iw\nBfA34GuqzM5Do2+0RwEnSyo1rR/4604Djq6UmbSg77Jp2kWn7yQcCZwqUttUpAin4EI2j1Ylk953\nuU5N9GFcltVrJJXNe/y82718BTAH+J8sNFUsp0XsnhUNOwNVnaPaY5KxUcA8VV2oql0LT4c1WqYR\nOvIuwM/Jngb0Ba4WYa1q1/oRwe3Af6hyVZ4aNdElOMf+I0nlG9XmeSUVkVTOwPUk9+tmE1AuOnOg\nI7SAGuno6QJ/GtoE4CcinF4tWEEEEeEM3KbFfbs5xSwTnZroNbiRS6eksnO1D0oqawF/xT27Tsk5\nGV1Hjvc2R95rBhsCz5T8vsi/FukBv2HoaGAxMFOEo0rnekVYS4Qf4nLDfEe1tqykTetK9DFgT1zq\n4mslle1K35dUtsad5HYKsIcmOqcIXZHaUeVRnA1PA64VYXkbCtviDl9yNlQeL0RXor/GRcFNklT+\n3T/4naZU+kgqR+Ay6j5FldFmpHH6dvemiEyCikcz/quq1nI6VUwh2wR+WP9VEQ7A7UP4lQgPA6uw\n7HSp0arML1RXok9KKmNwCfZul1RewTn6IcAA4OfAcTFZmF1UmSdClw3/IcIruI5HqQ2PaXY/Qd26\nEr1SUrkPtza1SFJ5GHgXl8JkIXC6JhoTHOaAaJOjLBG5A/iOqs6s8N4YYKKqTvC/fx9YqqpnV7g2\nOo5IJBJpAFVteG9SF92ODOqgmpAZwDAR2QR4FpeL5LhKF2bxx0QikUikMZoJLT1CRJ4BxgA3isjN\n/vXBInIjgKp+hDtg41ZcdMyVqvpY87IjkUgkkiVNTxNFIpFIpPUJvgPZ0qY0ERkqInf4zXQPi8gZ\n/vWBIjJJRJ4QkdtEpH/JZ77vtc8Rkf0K1NpHRB4QkesNa+wvIleJyGMi8qiIjDaq8/ve5rNF5HIR\nWdmCThG5SESWiMjsktfq1iUiI/3fNldEfl6Qzp95uz8oIleLlEQGGdJZ8t53RGSpyLLzma3pFJFv\n+O/0YRE5u+T1bHSqarB/QB9gHrAJsBIui+E2AfUMAnb2P68OPI47xOOnwHf962cCZ/mft/WaV/J/\nwzzgUwVp/TYu7fV1/neLGi8BTvU/98XtRDal05c1H1jZ/34lLqlecJ24oyWHA7NLXqtHV9fI/15g\nlP/5JmBCATr37fpecDvqTer0rw8FbgEWAAMt6gT2AiYBK/nf181aZ+iRgalNaar6vKrO8j+/BTyG\n2xdxKMtOdboEONz/fBhwhap+qKoLcYYYlbdOERkCHIjLF9O18G5N41rAeFW9CNz6kaq+bk0n8AYu\nE+eqItIXWBUX7BBcp6pOBl4te7keXaNFZANgDVW91193aclnctOpqpNUP0lKOB0XsmpOp+d/cRlU\nS7Gm82vAT/xzElV9MWudoZ2B2U1p4iKghuMq8vqqusS/tQQ+Oa94ME5zF0XpPweXpro0A6g1jZsC\nL4rIxSIyU0R+KyKrWdOpqq/gUho8jXMCr6nqJGs6S6hXV/nriym+jZ2K65lSQU9QnSJyGLBIVR8q\ne8uUTmAYsIeITBORTpFP8jhlpjO0MzC5ei0iq+O2vH9TVd8sfU/dmKs73bn+TSJyMPCCqj5AlZDe\n0Bo9fYERwAWqOgJ4G3fIyjIRBnSKyOa4jVeb4BrQ6iJy/HIiDOisWGjPuoIjIj8APlDVy0NrKUdE\nVgX+leWTL1oNce8LDFDVMbiO4J+zLiC0M1iMm6/rYijLe7PCEZGVcI7gMlW91r+8REQG+fc3AF7w\nr5frH+Jfy5PdgUNFZAEuWdfeInKZMY3g7LhIVe/zv1+Fcw7PG9O5CzBVVV9WFwp9NbCbQZ1d1GPn\nRSzbVVz6eiF6ReRk3HTmF0petqRzc1wn4EHfnoYA94vI+sZ04su+GsC3qaUisk6WOkM7g082pYlI\nP9ymtOtCiRERAX4HPKqq55a8dR3LTvo6CZe3pev1Y0Wkn4hsihvK3UuOqOq/qupQVd0Ud+j97ap6\ngiWNXufzwDMisqV/aR/gEeB6SzpxmS/HiMgq3v774PbEWNPZRV129nZ4Q1wklwAnlHwmN0RkAq4H\ne5jqcsecmtGpqrNVdX1V3dS3p0XACD8NZ0an51pgbwDfpvqp6kuZ6sxyFbzBlfMDcFE784DvB9Yy\nDjcPPwt4wP+bAAwE/g48AdwG9C/5zL967XOA/QvWuyfLoonMacQdnH4f8CCuV7OWUZ3fxTmq2bhF\n2ZUs6MSN/J4FPsCtrZ3SiC5gpP/b5gHnFaDzVGAuLqFcVzu6wJDO97u+z7L35+Ojiazp9HXyMl/u\n/UBH1jrjprNIJBKJBJ8mikQikYgBojOIRCKRSHQGkUgkEonOIBKJRCJEZxCJRCIRojOIRCKRCNEZ\nRCKRSIToDCKRSCQC/H9fP7MS26NF/AAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x7fcaf7495990>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "plt.plot(f_oscillation(1440))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Observe that we specified a parameter `a` above that we multiplied to the rotation matrix. Try setting it to `.999` and `1.001` respectively and observe the results."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We observe that oscillations increase to explosion or decrease quite rapidly to zero. Let us attempt to feed some input to modulate this behaviour.\n",
    "\n",
    "Write a function that takes a sequence of modulation factors. Multiply their exponential to the rotation matrix."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "# a_vector = T.fvector()\n",
    "# def controlled_rotation(\n",
    "\n",
    "# f_controlled_oscillation"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 160,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "# %load scan_tutorial_solutions/sol_theano_scan_modulated_oscillation.py\n",
    "a_vector = T.fvector()\n",
    "\n",
    "def controlled_rotation(a, old):\n",
    "    return T.exp(a) * rot.dot(old)\n",
    "\n",
    "controlled_oscillation, updates = theano.scan(controlled_rotation,\n",
    "                                             sequences=a_vector,\n",
    "                                            outputs_info=np.float32([1., 0.]))\n",
    "\n",
    "f_controlled_oscillation = theano.function([a_vector], controlled_oscillation)\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 158,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "[<matplotlib.lines.Line2D at 0x7fcaf6b59b50>,\n",
       " <matplotlib.lines.Line2D at 0x7fcaf6b59dd0>]"
      ]
     },
     "execution_count": 158,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYMAAAEACAYAAABRQBpkAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAIABJREFUeJztnXeYXFX5xz9fqoAKIv5oiYBAKIJKNdQs1UgvUkIRASV0\nbHTl5iICikqRIqD0KkVAaQmQACJSlBJKBKRIEQSkg5DA+/vjnIXNMrs7u/fOPXd238/z7JOdmTPn\nfLO7c997zttkZjiO4zhDmxlSC3Acx3HS48bAcRzHcWPgOI7juDFwHMdxcGPgOI7j4MbAcRzHoQRj\nIOkMSS9ImtzLmBMkPSrpPknLFV3TcRzHKZcydgZnAqN7elHSBsBiZrY4sBtwSglrOo7jOCVS2BiY\n2a3AK70M2QQ4O469A5hL0rxF13Ucx3HKowqfwYLA010ePwMMq2Bdx3Ecp0mqciCr22OvgeE4jlMj\nZqpgjWeB4V0eD4vPTYckNxCO4zgDwMy633APaJLCX8DCwOQeXtsAuCZ+PxL4aw/jDOyCMvS08gsY\nl1rDYNHZDhpdp+uc/jkbDvYi2LKp9XXRaWXMU3hnIOlCYBQwj6SngQyYOSo81cyukbSBpMeAt4Cd\ne5luFYn1zRhfVJfjOE6ZSAg4Efi1GT2G0rcrhY2BmY1pYszeTU63N3CyxBfNeLeYMsdxnFLZCBgB\nbJ1aSCuoVQayGVcD/yDkI9SVSakFNMmk1AKaYFJqAU0yKbWAJpmUWkCTTEotoEkmdX4jMSNwFLD/\nYL1RVTxzSo4kMzNJfAW4FljMjLdS63Icx5H4FrArsKZZvaIhO6+dheepmzEI33MRcL8ZRyaW5TjO\nEEdiZuBRYHszbkutpzuD3RgsCdwMLGLG22mVOY4zlJHYEdjJjHVTa2lEWcagVj6DTsyYAtwO7JRa\ni+M4Q5cYQbQ/cExqLa2mlsYg8gvg+9Fx4ziOk4KvESomDPpw9zobg9uAlwiF7hzHcVLwQ+AXdXMa\nt4LaGoP4w/81sEdqLY7jDD0kRgDLABen1lIFtTUGkcuAr0gsmlqI4zhDjt2As8x4L7WQKqhlNNH0\nz/MLYJoZByWQ5TjOEERiVkLp/VXNeCy1nt4Y1NFE3TgN+JbELKmFOI4zZNickOtUa0NQJrU3BmY8\nAjyMO5Idx6mO3Qg3okOG2huDyNnADqlFOI4z+JFYCPgScEVqLVXSLsbgcmAtic+mFuI4zqBnDHDp\nUHEcd9IWxsCM1wnF67ZKrcVxnEHP9sAFqUVUTVsYg8h5wI6pRTiOM3iR+BIwJ/Dn1FqqprAxkDRa\n0hRJj0o6sMHr80i6TtK9kh6Q9K0BLnU9MELiC4UEO47j9Mz2wAVmfJBaSNUUyjOQNCOhGc26hCb3\ndwFjzOzhLmPGAbOa2cGS5onj5zWzad3m6jNWVuIk4Dkzfjpg0Y7jOA2QmAF4Etiwndpa1iXPYGXg\nMTN70symAhcBm3Yb82/g0/H7TwMvdzcE/eBSYMsBvtdxHKc3VgVebydDUCZFjcGChCy9Tp6Jz3Xl\ndOCLkp4D7gP2K7DercAwPypyHKcFbAlcklpEKmYq+P5mzpgOAe41sw5JiwITJH3ZzN7oPjAeKXUy\nycwmTbeYMU3iCmALQolrx3GcwsS+BZsDG6fW0heSOoCOsuctagyeBYZ3eTycsDvoyqoQzvjN7J+S\nngCWAO7uPpmZjWtizUuBw3Fj4DhOeSwHTAUeSC2kL+JN8qTOx5KyMuYtekx0N7C4pIUlzQJsA1zV\nbcwUgoMZSfMSDMHjBdacCCwmTWeEHMdxirAFcPlQ6FvQE4WMQXQE700I+3wIuNjMHpY0VtLYOOxI\nYEVJ9wE3AAeY2X8HviZTgT8SfnmO4zhlsAWh0sGQpfYlrBuPZSPgh2bln5s5jjO0kFgKmAB8vh3z\nC+oSWpqKG4HlJT6TWojjOG3P5sAV7WgIyqQtjYEZ7xAcKKMTS3Ecp/3ZDPhDahGpaUtjEPkTsFFq\nEY7jtC8S8wIjCDlMQ5p2NgZXA6OlwuGxjuMMXUYDNwy1ctWNaFtjYMazhDoiqyaW4jhO+7Ih4cZy\nyNOW0UQfvYccmN2M/fu9Xi4RDMmGwKLAe4SciCsts9onnjiOUwyJmYH/AEuZ8XxqPQNlqEcTdfJH\nBuA3UK4VgL8AvwOmAVcSQss+B1yvXNco14gyhTqOUztWA/7ZzoagTNp9ZzADoSTGGmY81tR7co0F\nfgIcAJxjmX3Q7fVZgT0JNZV+YJmd0x9NjuO0BxLHAG+bUUo5h1SUtTNoa2MQ3sdvgclmHN/n2FwH\nAGOB9S2zf/YxdmlCxNKZltlP+qvLcZx6I/EgsLMZd6bWUgQ/JvqI64Cv9TVIuXYAdgfW7MsQAFhm\nDwGrANtGI+I4ziBBYmFgHhoUzByqDAZjcCOwusQnehqgXCsBxwKbWGbPNjuxZfYCsD6wezQmjuMM\nDjYErh3qWcddaXtjYMYrhLKzqzd6XblmB84D9hpIlFA0HpsCxyrXl4todRynNnwduCa1iDrR9sYg\nMp6ej4qOBP5mmf1+oJNbZpOBfYHLlOtTA53HcZz0SMwCrEGoouxEBosxuJ4GxiDeyY8B9im6gGV2\nIXAL8POiczmOk5RVgH+YMeBS+oORwWIM7iL0Rl6g84mYVHYskFtmL5e0zveBjZRrnZLmcxynetYn\n5BU5XRgUxsCMaQRH8vpdnt4EmBc4rbR1MnsV+A7wW+Warax5HceplPUIR8tOFwobA0mjJU2R9Kik\nA3sY0yHpHkkPSJpUdM0euJ5oDJRrBkKf5IMss2llLmKZXQf8DfpfAsNxnLRIfBZYErg9tZa6UcgY\nSJoROJFQ+W9pYIykpbqNmQs4CdjYzJYBvlFkzV64HlgvZiVvDLxPSBprBT8E9lOuz7dofsdxWsPa\nwC1epfTjFN0ZrAw8ZmZPmtlU4CJCGGZXtgMuM7NnAMzspYJrNsSMp4EXmWHq8sCPgCMsa016tWX2\nJMEIujPZcdoL9xf0QFFjsCDwdJfHz8TnurI4MLekiZLulrRjwTV743pW/M0ewBzAFS1cB4IhGKVc\ny7V4HcdxSkBCBH+BG4MGFG0M08yd98zA8sA6wOzA7ZL+amaPdh8oaVyXh5PMbFI/9dzAUn84C/hR\n9wJ0ZWOZvaVcRxGK3nnHNcepP4sTrnkPpxZSBEkdQEfZ8xY1Bs8Cw7s8Hk7YHXTlaeAlM3sHeEfS\nLcCXgY8ZAzMbV0jNjus9wbz3zcOde15WaJ7mORX4gXKtapn9paI1HccZGOsBE8yauomtLfEmeVLn\nY0mlVF0tekx0N7C4pIUlzQJsA1zVbcyVwOqSZpQ0O/BV4KGC6zZm0Rt25B+bPMc1J32lJfN3wzJ7\nlxC1dEQV6zmOU4j18ZDSHilkDMxsGrA3IZLnIeBiM3tY0lhJY+OYKYTKovcDdwCnm1npxkC5PgHs\nyr07/YFwJFUVZwOLKNcqFa7pOE4/iF3NRhHykZwGtH0/gw/fn2srYHfG2U+An5nx1fLU9bn2XoQe\nCd0jqRzHqQESqwEnmjHoAj68n8HH+SbhLv12YGmJuSpc+wzgq8q1TIVrOo7TPGsDN6UWUWcGhTFQ\nrnkJVQgvN+NdgkEYVdX6ltk7wPFAwwxsx3GSsxYwMbWIOjMojAGhMumVltmb8fFNhDuBKjkZ2EC5\nFqp4XcdxeiE2vlqZUHXY6YHBYgx2Aro2rr+Rap3IWGavAWcBe1a5ruM4fTISeNCM11MLqTNtbwyU\nawlCddKuW8C/AwtKzFexnJOAXZRrjorXdRynZ/yIqAna3hgAWwKXdc04NuN94GYqPiqyzB4H/gJs\nX+W6juP0ytq4MeiTwWAMvgFc2uD5G6nebwDBkbxvbK7jOE5CJGYHlgNuS62l7rS1MVCuLxAK4/25\nwcuV+w0inXcgKQyR4zjTsxpwrxlv9jlyiNPWxoBwRPQHy+z9Bq89DHxC4gtVCopls08A9qtyXcdx\nGuL+giZpd2PQ0xERsRjVTaTZHZwHrKZcw/sc6ThOK3Fj0CRtawyUaxiwKMFR3BM3Ev4YKsUyexu4\nENi16rUdxwlIfApYFm9x2RRtawyADYDrLbOpvYyZCKwdm1pUzenArspVtEy44zgDYw3gLjPeSS2k\nHWh3Y3B1bwPMeAJ4B1iqt3GtwDK7j9DvYXTVazuOA/gRUb9oS2OgXLMSftHXNzF8IgmOiiKnAbsl\nWttxhjpuDPpBWxoDQhG6Byyzl5sYm6JOUScXA6tH/4bjOBUh8RlgCeDO1FrahXY1Bn0eEXVhIjBK\nqv7/apm9BVwE7FL12o4zxFkT+GusYuw0QeELpKTRkqZIelRSjyWcJa0kaZqkLYquCWwIXNPMQDOe\nBV4GvlTCugPhNODbyjVjovUdZyiyFt6/oF8UMgaSZgROJDhJlwbGSPqYszaO+xmh/WWhyB7lGgHM\nBtzXj7clOyqyzO4F/pNqfccZori/oJ8U3RmsDDxmZk+a2VTCkUij1o/7EJLDXiy4HgTDc13M9G2W\nm0jnRIbQge1bCdd3nCGDxDzAwsDfEktpK4oagwWBp7s8fiY+9yGSFiQYiFPiU0WbLq8D3NDP90wC\n1pBIFfN/IbChcs2ZaH3HGUqMAm4zo7ccJKcbRS+OzVzYjwMOMjOTJHo5JpI0rsvDSWY2abrXQwLX\nKOA7/RJpvCjxL2B5EkQXWGYvKdeNwNaEZDTHcVrHoD4iktQBdJQ9b1Fj8CzQtf7OcMLuoCsrABcF\nO8A8wNclTTWzq7pPZmbj+lhvReApy+w/A9A6kXBunyrU7CzgINwYOE6rWRvYMbWIVhFvkid1PpaU\nlTFv0WOiu4HFJS0saRZgG2C6i7yZfcHMFjGzRQh+gz0aGYImWYdQb2ggpMw3gOA8XzQ6wB3HaQGx\nu+H8wL2ptbQbhYyBmU0D9iZkAj8EXGxmD0saK2lsGQK7sS4DNwY3A6tIzFKinqaJNZTOJ/Rrdhyn\nNXQAt8Ruh04/kPUrKKd1SDIz69mfkGt2Qojm/JbZGwNbg7uB75lx6wBlFkK5vgT8CVikhx4MjuMU\nQOJU4GEzjkutpSr6unY2SztlIK8O3DtQQxBJelRkmd1PCK9NGebq1AzlknLNo1zzecRZYToYxM7j\nVtJO5ZWL+As6mQgcCOTF5QyYswg5B/0Nj3UGCbE/9khCc6ZRwBcJ1XWnAnMo1zvA3wklVy62zF5I\npbWdkFiAEKQyObWWdqSddgYdFLf4fwZWlJituJwBcyGwkXJ9OqEGJwHKNbNy7QJMAc4EXgO+C8xj\nmc1tmc0LHzZkOZUQCj1Fuc5Urs+n0t1GdAA3m/FBaiHtSFv4DJRrDoK/YB7LrFCjConbgMPMCu8y\nBq4h15WE3s1npdLgVItyjSaUbnkKOBy4pZkseuWam2Aw9gR+Cfzc/U2NkTgdmGzGCam1VMlQ8xmM\nJPgLyuhY1JlvkJLzgB0Sa3AqQLnmUq4LgZOBfSyzdSyzm5stp2KZ/dcyO4yQr7M+cKNyzd9Cye1M\nB+4vGDDtYgzWgNIigFLnGwD8EVheuRbsc6TTtijXSoSz/5eAZSyzawc6l2X2FCG0ehJwu3ItXYrI\nQYLEMGAu4MHUWtqVdjEGa1KeMbgdWDY2y06CZfY/4HJgTCoNTmtRrp0IDuADLLN9LLO3i85pmb1v\nmY0DfgxMVK6Vi845iOjA/QWFqL0xUK5ZgJWA28qYLzbHvouw20iJHxUNQmKY6I+BccAoy+zSstew\nzM4FdgX+GHNXnGAMJiXW0NbU3hgQIioet8xeLXHO1CWtAW4BPqtcyybW4ZSEcs1AqM67ObCqZfZw\nq9ayzP5EyP6/Trm+0Kp12ohBXZyuCtrBGJTpL+gkuRPZMvuAUJ7CdweDgGgITiXkDIyyzP7d6jUt\ns0uAI4GrhnKossTnCSG5D6XW0s60izG4peQ57wQWl5i75Hn7y7nA9vFC4rQp8ff3G2BJ4OsFs+T7\ny0mE/JnzhvDfUQfBX1CPOPk2pdZ/PPGPezVK3hmY8R7wF4JjOhmW2YOE8hSjUupwBk7MJj6O0PZ1\nA8vszSrXjyGq+xAyb79b5do1ogM/IipMrY0B4QP2Sou23MmPiiLnMYhrrw8BDiBcjDaseEfwIbEi\n7nbAQcq1QgoNiVkLdx4Xpu7GYHXCFrgV1CHfAEJ5is2VK2WJjLZGQhKzSnyiibEzxbGzSD133Wtq\n3Vw7EjKDv26ZvVZkrqJYZk8C+wIXxoz9IYHEwsDsQMuc9UOFuhuDkYTjnFZwDzBMYt4Wzd8Ultlz\nhCZBG6fU0S5IzCixhsQ4ifESTxEKvL0BvCbxlsQ/Jf4Qx3xH4niJv0g8RygI9wbwdhz/gMQ5ErtL\nNJ0EqFzrAb8gGIJnW/Bf7TeW2UWEv6WUhRirpgOY5P6C4hQ2BpJGS5oi6VFJBzZ4fXtJ90m6X9Jt\nUr/iokcCdxTV2AgzphEc0x2tmL+fnItHFfWKxHwSPwGeAH4NfAI4nnBEMJsZs8Tn5gU2AC4gVOVd\njdCK9SBCvkrn2JmBhQhHdLfGcZMlbpbYQmLGHrWEbnXnA1tZZnWLYNkP2FG5VkwtpCI6cH9BKRQq\nVCdpRuAfhDT5ZwnJXGPMPoqvlrQK8JCZvSZpNDDOzEY2mGu6YkuxQNdTwGcss2kDFtmrfr4LLGnG\n7q2Yv2kduT4FPA0sbpm9mFJL3ZD4LHAIsDPhAn+aGfe3aK1ZgU2AHwKfIZQ7v6LrXWfsN/BX4FeW\nWS37Wcfjqx8CK0Z/wqAkHvM9CXzNjCmJ5SSjLoXqVgYeM7MnzWwqcBGwadcBZna72YfnqXcAw/ox\n992tMgSRWjiRo+PxGmDr1FrqQvQDbAc8AMwGLGPG3q0yBABmvGvGJYQd6V7AEcB4iUUAlGtGwo7g\nproagsh5hCq/e6QW0mIWBmYh3JA6BSlqDBYk3NF28kx8rid2JVz0muGrhDuwVjIZmDsWuUqNl6eI\nSMwJXEI42tnUjD3NeK6q9c0wMyYAXwHGA3dK7IDxE+CT1DyEM4abfhf4kXLNk1pPC1kL9xeURlFj\n0PQvQdJawC6ErXczjKTFxiAWtZpI+tIUABOALyjXYqmFpERiWYIT9AVgJTPuTKXFjKlmHAOsx5fP\nPpq359mbF5Yd0w5HLzGH5SIGtzO5A/cXlEZRY/AsMLzL4+GE3cF0RKfx6cAmZvZKT5NJGhe/ch5n\nNVrkPO5GXY6KOo/Ztk+tJRUSaxNam44zYy8z3k2tCYBxeovNdp6VSy6+j1PuP1dirtSSmmQcsNVg\nrH8V/QVDMr9AUkeXa+W40uYt6ECeiXBetw7wHKHMQ3cH8ucJMf07mFmPd/pdnSDKtQRwvWW28IDF\nNYnEksD1wMKpt5ux/v2FBEfykNr6SmxDiBLa2qw+H/CY/3E7cCrj7DRCtvGqwPpmvJxUXBMo177A\n+pbZRqm1lInEooRowGGpP7epqYUD2cymESonXk8oEnWxmT0saayksXHYYYTIjFMk3SOpmW1/y4+I\nuvAPghNqkYrW6427gfcJ/pIhQzQExwHr1skQRI4nJDT9xoz3CYldE4CJqXNUmuRUYFnlWiW1kJJx\nf0HJFM4zMLNrzWwJM1vMzI6Kz51qZqfG779tZp81s+XiVzMNOapwHgPBWUhNspHjbmBIOZIlNgdO\nIIQHtixSaCDEEM1RwG6dO7X493Iw8AfghrofGVlm7xJ6Lv80tZaS6WAIHhG1krpmIFe5M4B69Dfo\n5Dxgm9jUZ1Aj0UG4c92ghoZgaeBXwDe61xyKBmEc4e/mymbKYCTmbGCYcq2TWkgZRH/B2rjzuFRq\nZwxiXZUlgXsrXHYisHbRWjVlYJk9AUwhZNEOWiSWAC4Gxpjxt9R6uhL/Bi8htKyc3GhMNAjfA/4N\nnC/V77PUSczVyYAjY5XVdmdJ4D3gn6mFDCbq+Ae8AvBA7BNcFU8A7wJLVLhmb5wNfDO1iFYhMQ+h\nP/ChZtyYWk9X4sXyFOAuy+zM3sbG0OSdgP8jXGzrzMWEBjDrphZSAusAN7q/oFzqaAxWopqQ0g+J\nf1S1CDGNXAKsrVyfTS2kbGLNn98Dl5vx29R6GrAzodXqXs0MjuGv3wB2jv6PWhI76/2M4O9od9aF\net1EDAbqaAxWgCTHBrVwIgPEcsjXAGNSa2kB4+BDJ2ytUK4vEi6YW1tmbzX7PjNeALYATpX4Yqv0\nlcAFwKLK1bbRahIzEZz6N6XWMtioozFYkRBiWTUTgY4anf2eTTiCGDRIjCbceW8XwzRrg3LNTjhK\nOXAglUjNuJuQXX+xRC17U8TExmOooSHuB8sDz5rxfGohg426XPiADytCLgDVVyA04xngv0BdsjVv\nABaMUS1tT6z/dBbBELyQWE4jjicELfTqJ+iDs4D7gV+WIahFnAGMjLugdmQdwmfDKZlaGQOC1b+v\nxZVKe6NOR0XvE8JM296RHHdbZwEnmnFLYjkfQ7nGEI4e9iiS+R19T3sAo+vqP7DM3ibkdTRbI6xu\nrIP7C1pC3YzBCqQ5IuqkTvkGEI6Kdoilk9uZPYE5gKNTC+lOLAx4AsFPULiHsRmvEfoR/0ZigaLz\ntYhTgI2Va/7UQvpDPH4bCdycWstgpG7GYEXSOI87mQSsGZ1UyYmVJ5+nJruVgSAxguA03il2l6sN\nyjUrwU8wzjIrLa/FjL8CvyEYhNrF9VtmrxCKIiZt6jQAVgUmm/F6aiGDkboZg6Q7AzP+Q+jPsFwq\nDQ1oW0dyNKrnEKqQPpJaTwN+TuiUdXIL5v4poflKXSPCTgB2V666Z093xY+IWkjdjMF8pO9aVKd8\nAwh3cBsp16dTCxkA+wFv0ZqLbSGUazNCi8tdW1Eh1oz3CJFTx9axoJ1l9jBwD7Btai39wJ3HLaRu\nxuDe6DhNSW2cyACxJ/IkQmJT2yCxECGEcWzM1K0NyrUQoSbSGMvs1VatE8tsnAGc1Ko1CnIcsF87\nlKiIBQGXJpQTd1pA3YxBHWrU3AysIlGnQnFtVZ4inpOfBBxrxmOp9XQlFgC8CPiFZT331yiRHFhG\nYuMK1uov4wn9pddMLaQJOoDba9PwaBBSN2OQMpIIADNeAR4Fmim1XRVXA19Uri+kFtIkWxL6QxyT\nWkgDjiW01KwkF8CM/xFKW5wgMXsVazZLLFFxArBPai1N4P6CFuPGoDG1CjG1zN4Dzif0kK41sZn9\n8YTjofdS6+mKcu1EqGuzU7wQVkIsxncH9cz8PQ9Ypw3CTNfF/QUtpbAxkDRa0hRJj0pqmMgi6YT4\n+n2SeovUqUvESd2cyBB6SO+sXLUIe+2FccC1Zvw5tZCuKNdywC+AzWPtp6r5AbBHDLWtDZbZ64TC\niLW90Yj+p88SHN5OiyhkDCTNCJwIjCY4d8ZIWqrbmA2AxcxscWA3QsJLQ6q8W+uDW4GV6lRjJuYc\nPEWN+xzEftI7AIek1tKVWP31cmDPgdQdKgMzngWOBE6sYe7BqcBuNU5u/Bowvm6BCIONojuDlYHH\nzOxJM5tKcMxt2m3MJgQHKGZ2BzCXpNqF2nXFjDcINWbq1jf2NIJBrSu/Ao6K+Rq1IF7gzgcus8wu\nSSzn18D8wGaJdUyHZfY34D+Ei24d+Rqhz7rTQooagwUJSVqdPBOf62vMsILrVsFEgtOqTlwCrKpc\ntfv5SWwALErYKdaJccCswEGJdWDGVOD7wDESs6bW043fUMOMZImZCZ/D8am1DHaKGoNmk3W6b4vb\noUPReGp2pxRr7F9Ezc53Yxjur4Dv1clprFybErK3t0lY/HA6zJhAqMrbVPOcCrkIWF25hqcW0o2v\nAo/XrdKtch2hXHXepfebosbgWaDrH89wwp1/b2OGxec+hqRxXb46Cmoryl+AxWqYPXo6sGvNznf3\nAp4w45rUQjpRrsUJP6utLbPaHFtF9gcOju0/a0G80bgA+HZqLd2o6xHRSKY/8agMSR1dr5WlzWsF\nMvElzUQoH7EO8BxwJzDGzB7uMmYDYG8z20DSSOA4MxvZYC4zs1o51iQuA64045zUWrqiXHcDP7LM\nrkuuJVzQHgbWNOPhvsZXQay3cztwumVWu1IYABK/BjCrT4y/ci0LXAcsVJedlMRdwP5mTEqtpZOY\nsf0ysJRllnzHUta1s9DOwMymAXsTLPdDwMVm9rCksZLGxjHXAI9LeowQtbBnQc1Vci3w9dQiGnA6\n8J3UIiIHAZfUxRBEjiWEKfcYuVYDxgHbSizV18CqsMwmE3b266fWAh/eaIwg7NLrxMLA23UwBGVS\nOGbdzK4lXDS7Pndqt8d7F10nEdcBP5OYsWZtGi8EjlKuBS2zhkduVRC7l+0MLJNKQ3dio5p1gRVa\nUYCuLMx4WeIoQpb2Rqn1dOFMwu+0Dkd+6wGT6uSHiiwP/D21iLKpWwZyrYitMJ+hXqUpOhOFLgTG\nJpbyI+C3Zvw7sQ4AlGsEHzWqaYea9ycR6hatkVpIFy4C1ou5GakZTT39BSvgxmBIci31TPQ6kZAo\nlCREUWIxQiXVn6VYvzvKNRsh9PbHlllbZKrGomuHAUfVJREtVnG9htCtLRnx57E+YXdeN5anHkU1\nS8WNQd/U0m8Q69FPBrZKJCEHjjPjv4nW785xBL/VqX0NrBnnA3NSz6OilHwFeNOMxxPrmI7oPPZj\noiFKXUNMIWS0Vh6NIvElQgTZ8VWv3Qjl2o5QWHBsnf0EjYi+qEOAIyXqEi58EzCPcn0loYaNgT8m\nXL8nhgMfEKInBxVuDPogZo3eSM0S0CJXA/+nXFX7NH4CHB3LdiRFuZYgGKWt2sRP0Ig/Aa+T+Gim\nk9hg6izS7g42Ivxc6sZKwF3tdtPRDG4MmqOuR0XvE1pKVhatJTGS0CP6N1Wt2aOWXLMT/ASHWmb3\npdYzUMwwQoju4TUqU3EWsF1sBlQpEvMDixMKRtaNlahPqf1ScWPQHNcC68cG73Xjd8DGyjVfRev9\nFMhj05bUnEDwm5yeWkhRzLiV4PNIHSEGgGX2OPAgJOnQtgFwfdyV142VgLtSi2gFbgyaIJYffhJY\nPbGUj2GZ/ZdQRmDfVq8lsS7hzPTsVq/Vp5ZcOxJ+H7sPoi37IcAhEp9KLSSSypG8MTU8IlKuGQhh\npW4Mhjhm1eq2AAAZ+0lEQVRX8vHy3HXhl4Qw05ZdRGKo35HAj81IWqpAuZYiFMbb2jJL7rcoCzPu\nI/invp9aS+RSQvG6yoInJD5BCAa4tq+xCVgMeNUyezG1kFbgxqB5rgA2rUs8eFfilv4GWluiYlNg\nFsIZfTK6+AkOsszuT6mlRWTAvhLJk75i8bqrgG0rXHYt4H4zXq5wzWYZtEdE4MagP0wmlOKuTemF\nbhwDfK8VDr8Y8ngEcGgNuk2dSGh/eEZiHS3BjMeIxi61lsi5wI4VrrcR9QwphUHsPAY3Bk0TIz5q\ne1QUu1X9AxjTgunHAK+RuF5NbGi/CrDHIPITNOInwK7SxxpFpeAmYP54NNdS4q67riGl4DsDpwu1\nNQaRnwMHREdXKcTGNYcDB0eDmATlWprQ0H4ry+zNVDqqIAYs/I5Q+ymtlhC+fAHV7A6WA6ZCrSrg\nAqBcMwFfZhCWoejEjUH/uBX4Qk3u2BoxAXgb2LLEOXcFHjHjlhLn7BfKNQfh6OQAy+yBVDoq5mhg\nK4lFUwsBzgO2L/Mmowe2BC5LedPRC0sDz1hmr6UW0ircGPSDGEVzDbBJai2NiEcnhwF5GZ3QJGYn\n3J0eWnSuAWsItWBOIpzVnpVKR9VEB+oJhBpQabWEhL7XoOXVVbcALmvxGgNlUB8RgRuDgXAlsFlq\nEb1wHfAqsE0Jc+0F3G6WdGu8C7AisOcg9xM04lhgPYllUwuhxY5kiaWBOajvBXdQO4+hoDGQNLek\nCZIekTRe0lwNxgyXNFHSg5IekNTy5KgWcx0wUmLu1EIa0WV3MC6ecw4IiTkJvXp/XJa2fmsIhdKO\nJvgJ3kqlIxWx9tPRhEiu1FwAbBFLhbeCLYHLa3pEBL4z6JODgAlmNoKQLNMoHG4q8D0z+yKhifRe\nUusjE1qFGW8SzubrvDu4Efg3sH2BOX4AXJ2qnaVyzUnwE+wTy3UPVU4Blo81oZIRO+r9jdaV2t4S\nuLxFcxciGsClgHtTa2klRY3BJnxUmuBsGlwgzex5M7s3fv8mIVJggYLrpub3wNapRfRE3B38iOA7\n6PednMT/EY6IkpxXRz/BGcB4y+yiFBrqQqwBdTihJlRqzqMFR0XRST4fcFvZc5fECsBDltnbqYW0\nkqLGYF6zD5tCvwC91/yXtDAhfOyOguum5mpgldiwu5ZYZrcSGnB8bwBvPwi4wIwnSxXVPN8FPk99\nyjKk5ixgeKwNlZLLgTWV63Mlz7sFcEXN+ox3ZRVCX5NBTZ9nypImEKx2d6aLMDEzk9TjeZ+kTxJq\nnewXdwiNxozr8nCSmU3qS18KzHhL4jpgc+pdMXN/4A7lOsMye76ZN0gMB3YCvthSZT2tn2t1gjEa\naZm9m0JD3TBjqsRhhAY4N6Y6V7fM3lCuqwm74pNKnHorEkasNcEqJC7D0hVJHUBH6fNagQANSVOA\nDjN7XtL8wEQzW7LBuJkJWYXXmtlxPcxlZla7uj89IbElsLsZ66XW0hvKdQzwGcvs202NF6cBL5tx\ncGuVNVg71+eBvwK7WmZ1LFSWDIkZCDu93Iw/JNORazQwzjIrxYchfdi3YFjqAoiNiEeWzwGrWGZP\nJpbTkLKunUWPia4i3EUS/72i+wBJImRTPtSTIWhTrgVWiufrdeYIYCPlWqGvgfGDuTkhk7lSYgG6\nK4Bj3RB8nFgT6lDgiMTtMW8AFlauESXNtx1wcR0NQWSh+O9TSVVUQFFjcDSwnqRHgLXjYyQtIOnq\nOGY1YAdgLUn3xK/RBddNjhlvExLQysz2LZ2YMXkgcHoToaY/ITS5f6X1yj4i3n39jtBM5RdVrt1m\nXEPIIUnWHtMymwZcSPhMFyLWItoeOL/oXC1kFeD2oZDjUsgYmNl/zWxdMxthZuub2avx+efMbMP4\n/Z/NbAYz+4qZLRe/ritDfA24iJr0re2Dc4CX6cWZLLE8sCaQYvd2IKFW/G5D4UM3UKKv4GAgjzWj\nUnEusEM04kVYkXANqnP8/irA7alFVIFnIBfjWmCExGKphfRGvMCOBQ5Urp5q3RwJHGFGpcldyrUt\nIYx1c8vsnSrXbkdijahHgKZ8QC3iHuAdwq6/CNsB59c40QzcGDjNEHu0Xgh8M7WWvogNcI4Czu5+\nXCSxFuHOvNLIKOVai1B/ZwPL7Jkq125zDgV+FGtHVU68uShUniL6PbYlZDbXkpijszSDuFJpV9wY\nFOds4Jsx2qPuHEu4o/uwNHI8tz2K0M6ysgbkyrUscDGwjWU2uap1BwOxVtRtwD4JZZwPfEO5PjHA\n968DPGvGP0rUVDYrAg8OlR1rO1zA6s69wOvAqNRC+sIy+4CwixmrXJ0VKDcFPkG4MFdCjES5FtjP\nMptY1bqDjMOAH0p8rB5YFVhmTxP+9gdanmJX6t+tbsgcEYEbg8LE886z+CjEttZYZv8mVAK9SAf8\n30IEX8HBVbWzVK7FCLWTDrPMLqxizcFIrBn1R0INqVQM6KgoZu5/jRofEUVWZQhkHnfixqAczgc2\nk/hUaiHNEOP4jwO7hVlfe5lQibXlROf1TcDhllnd7wrbgRzYU+q9DEwLuQwYpVz9LcvyTeAqM15t\ngaZSiI181iAkxA0J3BiUgBkvABNpjzDTwDH/PonH1/0s+yxhjBt4qetmUa7lgFuAIyyzOpfwaBvM\neIpwI1J5tjiE8hSE3Ieme2dEH9W3gd+2SldJLA381zJ7LrWQqnBjUB4nE+7S2qOkxlvz/ZArz7iG\nT77wBnBukd4HfaFc6wDXA/taZqe1ap0hyk+BHSU+n2j9/h4VrQLMSP3vuEdBulavKXBjUB43Ehyx\nq6YW0hcSCwD7MW22AwgZ1HMD5xeIDGm8Ti4p1w8Id69bWWZ1bWnYtsRd6W+AcYkkTAAWUq4lmhw/\nFvhtzXMLICRg3pxaRJW4MSiJ6IA9BdgztZYm+ClwmhlPWmb/46OezjeUVZ5YuT5LqPS4LfBVy2xI\nfbAq5ufA1yX6rD9VNv0pTxF9G5sQSo/UlphZvSa+M3AKcDawQZ2L18WyE6MJuQUARIMwhnAndI9y\nfW3A84fdwFbAA8AzwBqW2aAv8pUSM14j5I6ckOiYsrM8RV/Xk7HA7834bwWairA4oUPjkPq7dWNQ\nIrHA26XA7qm1NCJeKH4FZGa83vU1y+wDy+xQQojs6cp1rnIt1GiehnMHI9C5tc6ALSyz70ZD47Se\nM4FZSBPEcC/wFr2Up4i1lHYnZJzXnTWBm4danSw3BuXzK2AviVY1Di/CNsBn6GWbbpndSGhs8zhh\nl3CRcn1dueboPjYagMWUa19CvZoz4txftsyGTLJOHYjHlPsCP5P4ZKVrf1SeoreyLFsDD5rxYDWq\nCjHkjoigYHObMmm35ja9IXEFMN6Mk1Nr6URiTuAhYCuz5hJplOszhBLD2xDalT5FOPr5APg0oZ7R\n+wQn4rnATTHL2UmExDnAM2YcUum6uYYB9wMLdN8Nxh3pPcChZlzd6P11QrmeAta3zOpcKuNDyrp2\nujFoARKrEpqHj6hL0w6J44HZzfjOgN4fdgaLAgsCAt4E/gk8N9S203UmRordD6xqxiOVrp3rBuBU\ny2y6FpESGwOHA8vXPYpIuRYmdNubv13+rsu6drY82WgoYsZfJJ4l9HZNXnIhOo23JSTSDAjL7C3C\nReb+snQ55WPGcxKHA6dJrF1VmZFIZ87Bh8Yg7gp+TCiP3g4X13WAie1iCMpkwD4DSXNLmiDpEUnj\nJfVYMEvSjLHD2R8Hul4bcgSQSWkNbiwVfAqh/tDLKbU4lXESIedl14rXvRxYs1t48nrAHJCub3M/\nWZdw7DnkKOJAPgiYYGYjCAlXB/Uydj/CefVQsrbjgecpUPO9JL5HKFt9VmIdTkWY8T7wHeDIeGxU\nzbqhPMXVhF1o565gHPDTincoAyKGxq5D6PM85ChiDDYhxNUT/92s0SBJw4ANCLVIBoVPoBnilvgQ\nYJzErCk0SCxFMNK7tMOH0SkPMyYTMpNPqjj34Bw+ugHaApidCsujF+RLwCuW2b9SC0lBEWMwr5m9\nEL9/AXqsnHgssD8MvYtRjNq5H9ij6rXj8dRZhKY1j1e9vlMLfkpw+n+rwjVvBIbpB8O+DBwN7B93\nKu3AugzRXQH04UCWNAGYr8FLh3Z9YGYm6WNHQJI2Av5jZvdI6uhLjKRxXR5OMrNJfb2nDTgQuFni\nwlhHpioOAt4g3B06QxAz/iexHTBJ4lYzHmv5mplNU66zeWmJXwH/NGur8/f1CP61WhOvpR2lzzvQ\n0FJJU4AOM3te0vzARDNbstuYIwlbxmkEh9angcvM7GPJKYMptLQ7EscAnzOr5g5NYhRwEbCiGc9W\nsaZTXyT2JeSLrF5Fa1Nt9q3VGHH1rZx+13L2ysL3tXq9MohFGl8Ehltmte2z0Iiyrp1Fjomu4qPu\nXjsBV3QfYGaHmNlwM1uE4FS6qZEhGAIcDqwjsUafIwsSi4GdD3zLDYET+TXwMuHYpvVcedb3+d9n\n/sV3F1m8kvXKYVXggXYzBGVSxBgcDawn6RFg7fgYSQtI6inLcChFE32IGW8Qmpef2cpuaBIzEwzB\n2WZc36p1nPYiBjPsQOjGt30r15LYDFiGTz6fU31oaxHWZ4iGlHbiGcgVIoWIKrPyPyQxYuQ3wHBg\nk7pkPjv1QWJZQtvR0Wb8rQXzDwPuBrZknP5OKF2yXDtE5yjXA8CultkdqbX0lzocEzn957vAKIkx\nLZj7B4QuUtu6IXAaEcNNxwJXSixa5twxeu0C4AQzbrPM3iFk3+9c5jqtQLkWAT4H3JVaS0rcGFSI\nGW8SOosdL7FSWfNK7EaoWLlR99LUjtMVMy4HfgJMkFiwjDnjrvQE4G2m90v8DthFuWYsY50WsiFw\n7VAvsujGoGLMuI+QHXqFRGEHm8R3CLVf1jaj9ttxJz1mnAqcCtwksXAJU+5P6GWwddfkRsvsHuAl\nQlZvndkI+FNqEanxQnUJMONKic8RPozrmtHvUrnxbiwnhO6uU0UMuTN4MONnEu8Af5bYMN6k9BuJ\nA4HdgFE97Ep/S7j5GT9wta1DuT5JNGSptaTGdwaJMOO3hDv6myXW6897Y1vNKwgRECOrLlXsDA7M\nOIFQu+oGie/0p2yFxKwSJxPCytc045kehl4ArKtc8xdX3BLWAe6wzIb88aobg4SYcRYh/+JsieMl\nPt3beImZJHYltBl8mHA3VmVWszPIMOMSYA1gL2C8xDJ9vSf26/groQTNKr3ls1hmrxFqEw2oj0YF\nbIofEQFuDJJjxiRCgaw5gMclfimxZmfrQonZJFaS+BHwCOFObDMzDjLj3WTCnUGDGVOAlQiJpBMl\nrpbYXmK4xIwSM0gMk9hOYjzwe+CXwDfMeK2JJU4CxirXzK37X/SfqGcTQuntIY8bgxpgxktmfBtY\nEfgf4YP2H4n3gFcJzc7nJZxrjjLjzmRinUGJGVPN+DWwEHAp8A1CqOU7wFTg74T2p+cAi5lxXrPN\naiyzycCjwOat0F6ADuCf7ZAHUQXuQK4RZjxJKAJ4aDy/nRV4t006RDmDADPeJtx8nAkfZrVbCbkr\nJwF7E3YVdeEbdOnKNtTxnUFNMcPM+J8bAiclccdQRhLjFcBiyvWlEuYqjHLNRNipXJZaS11wY+A4\nTsuxzKYSyqXsnVpLZE3gX5bZE6mF1AU3Bo7jVMVpwDe69UhOxRjapwNbJbgxcBynEiyzFwhn9El3\nB8o1G6EszPkpddQNNwaO41TJL4A9lGuOhBo2Be60zJ5LqKF2uDFwHKcyLLNHgVtI2+tgJ0KIrNMF\nNwaO41TNMcD3UyShKdcCwEgadGYc6gzYGEiaW9IESY9IGi9prh7GzSXpUkkPS3pI0siBy3Ucp92J\nDWSegJb09eiLnYDLLLO3E6xda4rsDA4CJpjZCODG+LgRxwPXmNlShLILDxdY03GcwUEOHFbl7iDm\nFuwOnFzVmu1EEWOwCXB2/P5sYLPuAyTNCaxhZmcAmNk0M2umlonjOIMYy2wS8BThTr0qNgKes8z+\nXuGabUMRYzCvmXVWzHyBUDunO4sAL0o6U9LfJZ0uafYCazqOM3j4MfBj5Zq1ovX2IpTFcBrQa20i\nSROA+Rq8dGjXB2ZmkhqVTZgJWB7Y28zuknQc4TjpsB7WG9fl4SQzm9SbPsdx2hfL7C/K9SChvPWJ\nrVxLuZYFlmEQ1CKS1EEoslfuvGYDK30jaQrQYWbPS5ofmGhmS3YbMx9wu5ktEh+vDhxkZhs1mM/M\nrOnmGo7jtD/K9RXgOmBJy+zVFq5zAXCvZfbzVq2RirKunUWOia7io/O+nWgQqmVmzwNPSxoRn1oX\neLDAmo7jDCIss3sJ15KsVWso12LAeoTaSE4PFDEGRwPrSXoEWDs+RtICkq7uMm4f4HxJ9xGiiY4s\nsKbjOIOPQ4EdlGupFs1/MHCyt7bsnQEfE5WNHxM5ztBFub4HbAisZ1l5FyXlWga4CVjCMnulrHnr\nRB2OiRzHccri18CclN8r+efATwerISgTNwaO4yTHMpsG7AwcqVwLlTGncm0ELA6cUsZ8gx03Bo7j\n1ALL7AFC/+9zi2YmK9enCZnGYy2z98rQN9hxY+A4Tp34GfBm/LcIvwLGW2Y3FZc0NHBj4DhObbDM\nPgB2ADZTrl0GMkd832rA98rUNtjxaCLHcWqHco0AJgH7WWZNZw0r1yjgUmBNy2xIFMX0aCLHcQYt\nltkjwNeBXyvX7s28R7k6COUmthkqhqBMfGfgOE5tidnDfwLuAH5gmb3UYMxMwH7AAcCYoeYnKOva\n6cbAcZxao1yfAg4HdiSUy78G+BfwaWB1YCzwb0Lk0GOpdKbCjYHjOEMK5VqE0Dt5FKGa8tvA34Fz\ngYllZi63E24MHMdxHHcgO47jOOXhxsBxHMdxY+A4juO4MXAcx3EoYAwkzS1pgqRHJI2XNFcP4w6W\n9KCkyZIukCprfu04juM0SZGdwUHABDMbAdwYH0+HpIUJ9cmXN7NlgRmBbQusmZzYjLr2tIPOdtAI\nrrNsXGc9KWIMNiEkgBD/3azBmNeBqcDskmYCZgeeLbBmHehILaBJOlILaIKO1AKapCO1gCbpSC2g\nSTpSC2iSjtQCqqSIMZjXzF6I378AzNt9gJn9l1Cf/F/Ac8CrZnZDgTUdx3GcFjBTby9KmkDI9OvO\noV0fmJlJ+lj2mqRFge8CCwOvAZdI2t7Mzh+wYsdxHKd0BpyBLGkK0GFmz0uaH5hoZkt2G7MNsJ6Z\nfTs+3hEYaWZ7NZivHqnQjuM4bUYZGci97gz64CpgJ0JHop2AKxqMmQL8WNJswP+AdYE7G03mpSgc\nx3HSUWRnMDfwe+DzwJPA1mb2qqQFgNPNbMM47gCCsfiAUFTq22Y2tQTtjuM4TknUplCd4ziOk47k\nGciSRkuaIulRSQcm1jJc0sSYJPeApH3j8z0m2MWkukfj/2H9ivXOKOkeSX+sq05Jc0m6VNLDkh6S\n9NW66WyUGFkHjZLOkPSCpMldnuu3LkkrxP/bo5KOr0jnMfF3fp+kyyXNWUedXV77gaQP4olHLXVK\n2if+TB+Q9LPSdZpZsi9CEtpjhGijmYF7gaUS6pkP+Er8/pPAP4ClgJ8DB8TnDwSOjt8vHTXPHP8P\njwEzVKj3+8D5wFXxce10EnJQdonfzwTMWSedcZ3HgVnj44sJx5rJNQJrAMsBk7s81x9dnTv/O4GV\n4/fXAKMr0Lle588FOLquOuPzw4HrgCeAueuoE1gLmADMHB9/rmydqXcGKwOPmdmTFvwIFwGbphJj\nZs+b2b3x+zeBh4EF6TnBblPgQjObamZPEn4RK1ehVdIwYAPgt0Cn871WOuPd4BpmdgaAmU0zs9dq\nprNRYuRzddBoZrcCr3R7uj+6vqoQ6fcpM+sM3DiHxgmipeo0swlm9kF8eAcwrI46I78itMzsSt10\n7gEcFa+TmNmLZetMbQwWBJ7u8viZ+FxyFEppLEf4Q+4pwW4BguZOqtR/LLA/wTHfSd10LgK8KOlM\nSX+XdLqkOeqk0xonRk6ok8Zu9FdX9+efpfrP2C6EO1Ma6EmqU9KmwDNmdn+3l2qlE1gcWFPSXyVN\nkrRi2TpTG4Naeq8lfRK4DNjPzN7o+pqFPVdvulv+f5K0EfAfM7uHj3YF04uogU7CsdDywMlmtjzw\nFt1qWKXWqekTIxcAPilph+kE1ONn+fFF+9aVHEmHAu+Z2QWptXRH0uzAIUDW9elEcvpiJuAzZjaS\ncBP4+7IXSG0MniWc13UynOmtWeVImplgCM41s87ciRckzRdfnx/4T3y+u/5hVFN7aVVgE0lPABcC\na0s6t4Y6nyHcdd0VH19KMA7P10jnisBfzOxlM5sGXA6sUjONXenP7/iZ+Pywbs9XolfStwhHmdt3\nebpOOhcl3ATcFz9Lw4C/SZq3ZjqJa18OED9PH0iap0ydqY3B3cDikhaWNAuwDSGZLQmSBPwOeMjM\njuvyUmeCHUyfYHcVsK2kWSQtQtjKNUyqKxMzO8TMhpvZIoQqsDeZ2Y411Pk88LSkEfGpdYEHgT/W\nSOcUYKSk2eLvf13goZpp7Eq/fsfxd/C6QhSXgB1pnCBaKpJGE+5gNzWz/3XTXwudZjbZzOY1s0Xi\nZ+kZQoXlF+qkM3IFsDZA/DzNYmYvlaqzTC/4AD3nXydE7TwGHJxYy+qEM/h7gXvi12hgbuAG4BFg\nPDBXl/ccErVPAb6WQPMoPoomqp1O4MvAXcB9hDubOeumk+A8fBCYTHDKzlwHjYRd33PAewTf2s4D\n0QWsEP9vjwEnVKBzF+BR4Kkun6OTa6Tz3c6fZ7fXHydGE9VNZ/ybPDeu+zdCKaBSdXrSmeM4jpP8\nmMhxHMepAW4MHMdxHDcGjuM4jhsDx3EcBzcGjuM4Dm4MHMdxHNwYOI7jOLgxcBzHcYD/B41at+aC\n/u3MAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x7fcaf6be4790>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "modulation = np.float32(np.concatenate([-np.ones(720), np.ones(720)]) / 500.)\n",
    "\n",
    "plt.plot(f_controlled_oscillation(modulation))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 2",
   "language": "python",
   "name": "python2"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 2
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython2",
   "version": "2.7.9"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 0
}
